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1  Introduction 

The  main  objective  of  the  proposed  research  has  been  to  investigate  the  proper  integration  of 
large-core  graded-index  multimode  fibers  (GIMFs)  in  ultrafast  mode-locked  fiber  lasers  to  dramat¬ 
ically  increase  the  pulse  energy.  Detailed  theoretical  investigations  are  aimed  equally  at  developing 
a  fundamental  understanding  of  pulse  propagation  in  GIMFs,  as  well  as  providing  key  opportuni¬ 
ties  to  enhance  the  performance  of  ultrafast  mode-locked  fiber  lasers. 

2  Research  Outline 

In  the  following,  the  main  results  of  the  research  are  outlined  in  chronological  and  logical 
order.  Specifically,  the  rationale  behind  the  research  direction  as  captured  in  the  sequence  of  the 
publications  is  explained. 

2.1  Preliminary  analysis  of  nonlinear  light  propagation  in  GIMFs: 

In  the  first  publication  [ J 1  ] ,  a  detailed  analysis  of  the  modal  properties,  dispersive  behavior,  and 
nonlinear  mode  coupling  in  GIMFs  were  presented.  A  simplified  form  of  a  generalized  nonlinear 
Schrodinger  equation  was  derived  to  explore  the  rich  nonlinear  dynamics  related  to  the  propaga¬ 
tion  and  interaction  of  light  pulses  in  graded  index  multimode  fibers  in  a  tractable  manner.  Issues 
related  to  four-wave-mixing  (FWM)  in  these  fibers  were  briefly  discussed,  but  detailed  investiga¬ 
tions  on  the  FWM  properties  of  GIMFs  were  left  for  future  publications  (see  [J8]).  Publication 
[Jl]  was  mainly  concerned  with  laying  out  the  basic  formalism  and  was  meant  to  pave  the  way  for 
future  studies  on  detailed  treatment  of  specific  nonlinear  phenomena  such  as  self-phase  modulation 
(SPM),  cross-phase  modulation  (XPM),  and  FWM  in  GIMFs. 

2.2  Nonlinear  switching  based  on  multimodal  interference: 

The  results  obtained  in  [Jl]  and  earlier  publications  co-authored  by  the  PI-Mafi  showed  that  a 
multimode  interferometric  device  based  on  single-mode  fiber  (SMF)  and  a  multimode  fiber  (MMF) 
junctions  in  the  form  of  a  SMF-GIMF-SMF  coupler  can  be  used  as  a  nonlinear  switch  to  mode-lock 
fiber  lasers.  A  similar  device  geometry  had  previously  been  used  by  PI-Mafi  and  collaborators  to 
obtain  ultra-low-loss  coupling  between  two  SMFs  with  very  different  mode-field  diameters.  In  the 
second  publication  [J2],  the  authors  presented  a  detailed  investigation  of  the  nonlinear  multimodal 
interference  in  the  SMF-GIMF-SMF  coupler.  The  results  indicated  the  potential  application  of  this 
simple  geometry  for  nonlinear  devices,  such  as  in  nonlinear  switching,  optical  signal  processing, 
or  as  saturable  absorbers  in  mode-locked  fiber  lasers.  Saturable  absorption  in  this  all-fiber  config¬ 
uration  was  shown  to  provide  attractive  properties  to  be  used  in  high  pulse  energy  mode-locked 
fiber  lasers. 

2.3  Why  is  the  switching  power  so  high?  The  case  for  multicore  fibers: 

The  nonlinear  switching  behavior  presented  in  [J2]  provided  a  very  interesting  possibility  for 
building  an  in-line  all-fiber  nonlinear  switch.  However,  it  was  observed  that  a  much  higher  power 
was  required  for  nonlinear- switching  in  multimode  fibers  compared  to  the  previously  reported 
nonlinear  switching  based  on  the  nonlinear  mode-coupling  in  multicore  fibers.  Moreover,  the 
nonlinear-switching  quality  of  the  multimode  fiber  was  shown  to  be  inferior  to  those  reported  for 
multicore  fibers.  This  prompted  the  team  to  investigate  the  main  differences  in  nonlinear  switching 
between  multicore  versus  multimode  waveguide  couplers  in  publication  [J3]. 
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In  publication  [J3],  it  was  shown  that  the  switching  power  threshold  in  nonlinear  multimode 
junctions  was  larger  than  in  multicore  junctions,  mainly  because  the  value  of  the  modal  propaga¬ 
tion  difference  (supermodes  in  multicore  fibers)  was  substantially  larger  in  the  former  case.  The 
XPM/FWM  terms  also  played  a  role  in  setting  the  power  threshold  but  their  importance  was  or¬ 
ders  of  magnitude  lower  than  that  of  the  separation  between  effective  propagation  constants  of  the 
modes.  It  was  shown  that  the  injection  of  the  right  combination  of  the  modes  at  the  input  and  the 
collection  of  the  right  combination  of  the  modes  at  the  output  were  the  main  factors  behind  the 
quality  of  the  nonlinear  switching  behavior.  While  creating  the  proper  injection  profile  is  often 
trivial  in  a  weakly-coupled  multicore  waveguide,  obtaining  the  right  combination  of  the  modes  in 
a  multimode  junction  is  difficult,  if  not  impossible,  to  achieve.  Therefore,  multimode  junctions 
show  lower  quality  nonlinear  switching  behavior. 

2.4  How  many  cores  are  optimal  for  nonlinear  switching  in  multicore  fibers? 

Based  on  the  results  obtained  in  publication  [J3],  the  team  decided  that  it  was  important  to 
conduct  a  detailed  investigation  of  the  nonlinear  switching  behavior  of  multicore  fibers.  Several 
designs  had  been  suggested  and  even  fabricated  by  other  research  groups,  but  none  had  resulted  in 
a  robust  high-quality  nonlinear  switching  behavior.  At  this  point,  the  main  question  for  the  team 
was  whether  increasing  the  number  of  the  cores  improved  the  saturable  absorption  characteristics 
of  these  nonlinear  multi-core  arrays.  The  team  studied  the  saturable  absorption  characteristics  of 
two-,  three-,  and  five-core  one-dimensional  fiber  coupler  arrays  and  the  seven-core  hexagonal  fiber 
coupler  array  in  publication  [J4] .  It  was  shown  that  the  performance  of  all  these  saturable  absorbers 
were  comparable  and  not  much  was  gained,  if  anything,  by  going  from  a  two-core  nonlinear  cou¬ 
pler  geometry  to  a  higher  number  of  cores. 

2.5  Concentric  ring  core  fiber  for  nonlinear  switching: 

Previously  in  publication  [J4],  the  team  showed  that  going  from  a  two-core  fiber  geometry  to 
a  higher  number  of  cores  does  not  improve  the  nonlinear  switching  performance  of  the  device 
considerably,  if  any  at  all  [7].  Hence,  a  two-core  fiber  coupler  with  optimized  parameters  seemed 
to  be  the  optimum  solution  for  a  nonlinear  switch.  In  publication  [J6],  the  team  proposed  a  novel 
concentric  ring  core  fiber  geometry  as  a  nonlinear  switch  to  initiate  the  mode-locking  of  the  fiber 
lasers.  This  fiber  has  two  concentric  cores  composed  of  a  circular  core  located  at  the  center  of  the 
fiber  and  a  ring  core  placed  around  the  central  core.  The  main  reason  for  choosing  this  geometry 
is  its  simple  fabrication  process.  Given  that  conventional  fibers  are  fabricated  radially,  maintaining 
the  sensitive  specifications  of  this  fiber  during  the  drawing  process  is  easier  than  non-concentric 
conventional  two-core  fibers.  Another  related  reason  is  that  aligning  the  central  core  of  this  fiber 
with  the  SMFs  in  the  laser  cavity  is  much  easier  and  more  robust  compared  with  conventional 
two-core  geometries,  especially  in  the  presence  of  the  inevitable  variations  in  the  fiber  geometry 
due  to  the  fabrication  and  drawing  process. 

2.6  Multicore  design  guidelines 

In  publication  [J7],  the  team  has  reported  important  design  rules  regarding  multicore  nonlinear 
switches.  The  design  rules  will  set  limits  on  the  maximum  fabrication  errors  that  are  allowed  for 
the  nonlinear  switching  to  be  observed  in  multicore  fibers. 
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2.7  Extreme  Raman  cascade  generation  from  523  to  1750  nm  using  a  GIMF 

Four-wave  mixing  in  GIMFs  is  another  nonlinear  effect  that  has  been  of  interest  to  the  team, 
especially  for  the  generation  of  new  wavelengths  based  on  multimode  phase  matching.  In  this  line, 
it  was  decided  to  collect  some  experimental  data  on  FWM  in  GIMFs  before  proceeding  to  explore 
the  theory.  The  experimental  work  was  conducted  by  a  graduate  student,  Hamed  Pourbeyram.  The 
observations  were  totally  surprising.  The  team  reported  in  [J5]  on  the  generation  of  a  Raman  cas¬ 
cade  spanning  the  wavelength  range  of  523  to  1750  nm  wavelength  range  (perhaps  even  higher  but 
not  measured),  in  a  standard  telecommunication  graded-index  multimode  optical  fiber.  Despite  the 
highly  multimode  nature  of  the  pump,  the  Raman  peaks  were  generated  in  specific  modes  of  the 
fiber,  providing  substantial  beam  cleanup  during  the  stimulated  Raman  scattering  process.  The  ob¬ 
served  results  and  the  theoretical  investigation  that  have  followed  [J5]  clearly  establish  the  Raman 
versus  FWM  regimes  and  what  the  conditions  need  to  be  in  order  to  generate  FWM  without  much 
Raman  contamination  in  GIMFs.  Moreover,  the  team  has  been  looking  at  the  generated  Raman  for 
various  device  applications,  including  an  inexpensive  tunable  fiber-based  light  source.  The  green 
pump  laser  used  in  the  experiment  is  available  in  most  optics  laboratories  and  an  inexpensive  com¬ 
mercial  fiber  coupled  with  a  few  other  optical  elements  can  potentially  result  in  a  very  inexpensive 
high-power  tunable  light  source. 

3  Upcoming  results 

3.1  Four- wave  mixing  in  graded-index  multimode  fibers 

In  publication  [J8],  general  guidelines  are  established  for  the  generation  and  observation  of 
FWM  in  GIMFs.  Once  the  theoretical  work  is  finalized,  experimental  results  are  expected  to 
follow. 

3.2  Experimental  verification 

The  results  reported  so  far,  especially  on  nonlinear  switching  in  multimode  and  multicore 
fibers,  have  been  followed  with  interest  in  the  research  community  and  experimental  groups  have 
contacted  the  team  to  fabricate  and  test  some  of  the  designs.  We  expect  to  report  more  on  the 
experimental  verification  of  the  theoretical  results  in  the  future. 

4  Key  Personnel 

The  main  project  participants  are: 

•  Dr.  Arash  Mafi,  Associate  Professor,  Department  of  Physics  &  Astronomy,  Optical  Science 
and  Engineering,  Center  for  High  Technology  Materials,  University  of  New  Mexico 

•  Dr.  Elham  Nazemosadat,  Postdoctoral  Scholar,  Electrical  Engineering  Department,  Univer¬ 
sity  of  Wisconsin-Milwaukee. 
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Pulse  Propagation  in  a  Short  Nonlinear  Graded-Index 

Multimode  Optical  Fiber 

Arash  Mafi,  Member,  IEEE,  Member,  OSA 


Abstract — We  present  a  detailed  analysis  of  the  modal  proper¬ 
ties,  dispersive  behavior,  and  nonlinear  mode  coupling  in  graded- 
index  multimode  fibers  (GIMFs),  and  lay  out  a  simplified  form  of 
a  generalized  nonlinear  Schriidinger  equation,  which  can  be  used 
to  explore  the  rich  nonlinear  dynamics  related  to  the  propagation 
and  interaction  of  light  pulses  in  GIMFs  in  a  tractable  manner. 
We  also  briefly  discuss  an  application  of  the  presented  formalism 
in  the  study  of  four-wave  mixing  in  these  fibers.  While  the  reported 
formalism  is  fairly  general,  our  presentation  is  mainly  targeted  at 
device  applications  in  which  short  segments  of  GIMFs  are  used. 

Index  Terms — Fiber  nonlinear  optics,  four-wave  mixing,  optical 
fibers. 


I.  Introduction 

THE  traditional  view  of  multimode  optical  fibers  (MMFs) 
as  low-quality  waveguides  is  gradually  changing,  as  some 
researchers  over  the  past  few  years  have  begun  to  view  them  as 
the  new  landscape  of  opportunities  in  novel  applications.  For  ex¬ 
ample,  there  has  been  a  surge  of  interest  recently  in  multimode- 
multiplexed  fiber-optic  communications  for  ultrahigh-data-rate 
transmission  [1].  Recently,  a  pioneering  work  on  the  integration 
of  MMFs  in  mode-locked  fiber  lasers  was  reported;  MMFs  were 
used  to  reduce  the  total  nonlinear  phase  shift  in  the  fiber  laser 
cavity  in  order  to  increase  the  pulse  energy  [2].  We  also  recently 
reported  using  graded-index  multimode  fibers  (GIMFs)  to  build 
mode-field  adapters  and  ultralow-loss  couplers  using  the  multi- 
mode  interference  (MMI)  effect  [3],  [4]. 

While  nonlinear  phenomena  in  MMFs  have  been  studied 
since  the  beginnings  of  the  field  of  nonlinear  optics  [5]— [8], 
a  systematic  and  comprehensive  analysis  of  the  nonlinear 
behavior  of  light  propagation  in  multimode  has  been  lacking. 
Perhaps,  this  can  be  partially  blamed  on  the  complexity  of  the 
nonlinear  light  propagation  in  MMFs.  Even  with  only  a  few  key 
dispersive  and  nonlinear  processes  governing  the  propagation 
of  light  in  single-mode  fibers  (SMFs),  their  interplay  results  in 
such  vast  and  rich  dynamics  that  has  been  the  subject  of  intense 
research  over  the  past  few  decades  where  new  phenomena  are 
still  being  discovered  [9];  perhaps  there  has  been  little  incentive 
to  add  yet  another  complication  to  the  system  and  cross  over  to 
the  multimode  side. 
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Conventional  multimode  fibers  are  often  prohibitively  com¬ 
plex  to  analyze  and  their  rich  nonlinear  dynamics  are  buried  un¬ 
derneath  their  hundreds  of  propagating  modes  and  millions  of 
nonlinear  coupling  terms  [10].  In  this  paper,  we  will  argue  that 
GIMFs  have  specially  interesting  properties  that  make  them  an 
ideal  platform  to  explore  the  rich  nonlinear  dynamics  related 
to  the  propagation  and  interaction  of  multiple  modes  in  optical 
fibers  in  a  more  tractable  manner.  In  particular,  we  will  show 
that  the  nonlinear  dynamics  in  GIMFs  can  be  considerably  sim¬ 
plified  for  better  fundamental  understanding  and  easier  control 
in  novel  device  applications.  Our  arguments  are  in  line  with  our 
recent  work  in  the  “linear”  regime  on  making  an  ultralow-loss 
coupler  based  on  MMI  in  GIMFs  [3].  We  will  show  that,  with 
sufficient  care,  nonzero-angular-momentum  modes  are  never 
excited;  therefore,  the  nonlinear  propagation  dynamics  can  be 
limited  to  only  a  handful  of  zero-angular-momentum  modes. 
Moreover,  unlike  conventional  multimode  fibers,  all  the  modes 
in  GIMFs  can  propagate  with  nearly  identical  group  velocities  at 
special  wavelengths;  therefore,  nonlinear  coupling  among  short 
pulses  is  maximally  achieved  all  along  the  fiber,  adding  to  the 
richness  of  the  nonlinear  dynamics. 

The  presence  of  multiple,  though  small  number  of  nonlinearly 
coupled  modes  gives  rise  to  an  abundance  of  novel  and  complex 
nonlinear  phenomena,  above  and  beyond  that  of  the  highly  inter¬ 
esting  landscape  of  two-polarization-mode  fibers  [5].  The  inten¬ 
tion  of  this  paper  is  to  lay  out  a  comprehensive  formalism  for  the 
study  of  nonlinear  light  propagation  in  GIMFs  in  the  scalar  ap¬ 
proximation  (polarization  effects  are  ignored).  In  particular,  the 
modal  properties,  dispersive  behavior,  and  nonlinear  coupling 
terms  are  derived  and  explored  in  sufficient  detail,  paving  the 
way  for  future  detailed  studies  of  nonlinear  phenomena  such  as 
self-phase  modulation  (SPM),  cross-phase  modulation  (XPM), 
and  four- wave  mixing  (FWM)  in  GIMFs.  While  the  reported 
formalism  is  fairly  general,  our  presentation  is  mainly  targeted 
at  device  applications  in  which  short  segments  of  GIMFs  are 
used,  such  as  those  reported  in  [3]  and  [4]. 

Section  II  is  primarily  a  review  of  the  modal  and  dispersion 
properties  of  GIMFs.  While  most  of  the  material  presented  in 
this  section  can  be  found  in  various  earlier  publications,  our  in¬ 
tention  has  been  to  provide  a  relatively  self-contained  manu¬ 
script  with  a  more  modern  and  accessible  notation.  We  note  that 
the  modal  and  dispersion  properties  reviewed  in  Section  II  are 
the  key  ingredients  which  will  be  used  to  construct  the  equations 
governing  the  nonlinear  propagation  of  light  in  GIMFs;  there¬ 
fore,  a  more  detailed  treatment  of  these  parameters  is  warranted. 

In  Sections  III  and  IV,  we  construct  a  simplified  generalized 
nonlinear  Schrodinger  equation  (GNLSE)  for  nonlinear  light 
propagation  in  GIMFs.  In  Sections  V  and  VI,  we  apply  the 
GNLSE  in  a  couple  of  specific  examples,  such  as  FWM  in 
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GIMFs,  and  conclude  in  Section  VII.  Some  issues  related  to 
mode  coupling  and  excitation  of  nonzero-angular-momentum 
modes  are  briefly  discussed  in  the  Appendix. 

II.  Brief  Overview  of  GIMFs 

In  this  section,  we  present  a  brief  overview  of  GIMFs  and 
establish  the  notation  that  will  be  used  in  the  rest  of  this  paper 
[11].  The  refractive  index  profile  of  a  GIMF  is  given  by 

n2(p)  =  «o  [i  -  2A  (^)  ];  P<R 

=  n20(  1-2A),  p>R  (1) 


depend  on  the  Ge  doping  concentration.  In  particular,  the  Sell- 
meier  coefficients  are  available  in  [13]  and  [14]  for  pure  fused 
silica  (refractive  index  represented  as  ns;( A)),  and  Ge-doped 
fused  silica  with  a  Ge  doping  concentration  of  6.3-wt%  (refrac¬ 
tive  index  represented  as  ng. 3(A)).  Here,  we  have  implicitly  as¬ 
sumed  that  the  refractive  index  of  the  Ge-doped  fused  silica  in¬ 
creases  linearly  with  the  Ge  doping  concentration  and  we  have 
already  relied  on  this  assumption  to  relate  (1)  and  (7).  This  is 
a  reasonable  approximation,  which  is  also  supported  by  experi¬ 
mental  measurements  [15].  Therefore,  we  derive  a  generalized 
Sellmeier  equation  to  incorporate  the  dependence  of  the  refrac¬ 
tive  index  on  both  the  wavelength  and  Ge  doping  concentration 
of  fused  silica  as 


where  R  is  the  core  radius,  no  is  the  maximum  refractive  index 
in  the  center  of  the  core,  A  is  the  relative  index  step,  and  a  « 
2  characterizes  a  near  parabolic-index  profile.  The  transverse 
electric  field  of  a  propagating  mode  can  be  decomposed  into 
separate  radial  and  angular  components  as 

Ep,m{p,<j>)  =  Ep,m(p)eim*  (2) 

where  [12] 

and 

P0  =  (4JV2)V4>  Ap  =  \] 7r(p  -f'|m|)! '  (4) 


l\',"  is  the  generalized  Laguerre  polynomial.  Here,  each  mode 
is  labeled  by  two  integers,  p  and  m,  referred  to  as  the  radial 
and  angular  numbers,  respectively.  The  coefficients  in  the  field 
profile  in  (3)  are  chosen  such  that  the  fields  are  orthonormal 
according  to 

/‘2ir  /‘Oc 

/  d(f>  pdpE*,rn,(p.4>)EPtrn{p,(j))  =  8  ppf  brrim'  • 

J  0  Jo 

(5) 


n2{\,  d)  =  n2s j(A)  +  ^  [126.3(A)  -  n^A)]  (8) 

where  d  is  again  the  Ge  doping  concentration  in  wt%.  The  utility 
of  (8)  is  that  it  provides,  along  with  (1)  and  (7),  a  reliable  esti¬ 
mate  of  the  refractive  index  for  a  broad  wavelength  span  from 
near-ultraviolet  through  mid-infrared  region  for  common  levels 
of  Ge  doping  concentration  in  GIMFs. 

All  the  modes  with  equal  group  mode  number  g  =  2p+\m\  + 
1  are  almost  degenerate  in  the  value  of  the  propagation  constant, 
which  is  given  according  to  the  formula: 

0g  =  n0k(l  -  2AX)1/2  (9) 


where 


(  a  \2q/(q+2) 

-bfc)  ■  ,10) 

The  group  velocity  of  the  modes  in  mode  group  number  g  is 
given  by 


Vf  =  (P^1  =  Wg)-1  (ID 

where  (3g  defined  in  (9).  After  a  few  algebraic  steps,  one  can 
show  that 


The  total  number  of  guided  modes  (counting  the  polarization 
degeneracy)  is  given  as 


(i2) 

•y  c  \  a  +  2  J 


Na  =  ~^—n2k2R2  A  (6) 

q  +  2 

where  N2  in  (3)  is  also  given  by  (6),  and  k  =  2n/X. 

The  near  parabolic  shape  of  the  refractive  index  profile  in  (1) 
is  due  to  the  near  parabolic  Germanium  (Ge)  doping  concentra¬ 
tion  d  in  the  core  of  the  GIMF 


where  rt®r  =  dk(n0k)  is  the  group  index  at  the  center  of  the 
core,  c  is  the  speed  of  light  in  vacuum,  and  y  is  defined  as  [13], 
[16] 


y  = 


n0  k  d  A 
77 q1  A  dk 


(13) 


d{p)  =  d0 


—  0,  p  >  R 


P  <  R 


(7) 


where  d0  is  the  maximum  doping  in  weight  percentage  (wt%) 
at  the  center  of  the  core.  At  each  given  value  of  Ge  doping  con¬ 
centration,  the  refractive  index  can  be  accurately  described  by 
a  fit  to  the  Sellmeier  equation,  where  the  Sellmeier  coefficients 


There  are  three  important  points  of  caution  regarding  the  use 
of  analytical  equations  presented  previously.  The  first  point  is 
that  (3)  is  only  strictly  valid  for  a  —  2;  however,  the  deviations 
are  small  and  (3)  is  often  very  reliable.  The  second  point  is  that 
(3)  is  only  a  solution  of  Maxwell’s  equations  in  the  core  of  the 
fiber  and  must  be  matched  to  an  exponentially  decaying  modi¬ 
fied  Bessel  function  of  the  second  kind  in  the  cladding.  In  calcu¬ 
lating  the  overlap  integrals  for  coupling  in  and  out  of  the  GIMF, 
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wavelength  (|im) 

Fig.  1.  Relative  delays  rsi  =  Pg1^  —  Pi1^  (delay  relative  to  mode  group  g  —  1) 
are  plotted  as  a  function  of  the  wavelength  for  different  mode  groups. 
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Fig.  3.  Same  as  Fig.  2  except  at  A  =  1000  nm. 
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Fig.  2.  Relative  delay  is  plotted  at  A  =  850  nm  for  all  the  modes  in  a  GIMF. 
The  horizontal  axis  represents  the  mode  group  number  g  and  the  dots  represent 
the  values  of  relative  delay  calculated  using  the  FEM,  where  the  cladding  effect 
is  now  taken  into  account.  The  large  black  dots  are  for  m  =  0  modes,  the  small 
red  dots  belong  to  m  ^  0  modes,  and  the  green  dashed  line  represents  the 
analytical  expression  from  (12). 


as  long  as  the  beam  diameters  of  the  in-coupling  and  out-cou¬ 
pling  “often”-Gaussian  fields  are  smaller  than  the  core  diameter 
of  the  GIMF,  it  is  perfectly  fine  to  use  (3).  For  larger  in-cou¬ 
pling  beam  diameters,  the  use  of  (3)  can  still  be  at  least  quali¬ 
tatively  valid.  In  such  cases,  better  accuracy  can  be  achieved  if 
the  cladding  profile  is  properly  accounted  for,  for  example  by 
using  a  numerical  solver.  In  either  case,  small  adjustments  in 
the  normalization  coefficients  of  (3)  are  expected  according  to 
(5)  in  order  to  compensate  for  substituting  (3)  for  the  modified 
Bessel  function  of  the  second  kind  in  the  cladding.  The  third 
point  is  that  the  presence  of  the  cladding  removes  the  degen¬ 
eracy  among  modes  within  a  mode  group,  both  in  the  propaga¬ 
tion  constant  and  the  group  velocity.  While  this  effect  is  gener¬ 
ally  quite  small,  it  is  more  pronounced  for  higher  order  mode 
groups,  which  have  a  larger  overlap  with  the  cladding.  Numer¬ 
ical  examples  of  this  last  point  are  explored  later  in  Figs.  2  and 
3. 

A.  Optimum  Value  of  a 

The  optimum  value  of  a  is  often  chosen  to  minimize  the  mul¬ 
timode  dispersion  in  an  optical  fiber  communication  link  at  a 
specific  wavelength.  Olshansky  and  Keck  [16]  have  derived  a 


compact  formula  for  the  pulse  dispersion  in  GIMFs  and  have 
shown  that  the  minimum  pulse  dispersion  occurs  for 


Q^opt 


=  2  +  t/- A 


(4  +  y)  (3  +  y) 
(5  +  2y ) 


(14) 


While  the  optimum  value  of  a  from  (14)  is  a  classic  result 
and  is  often  quoted  as  the  optimal  design  point  for  high-band- 
width  GIMFs,  it  may  not  be  suitable  for  most  modern  multi- 
mode  optical  communication  systems.  The  key  assumption  in 
the  derivation  of  aopt  by  [16]  is  the  the  equipartition  of  the 
power  among  all  the  guided  modes.  However,  this  assumption 
heavily  favors  the  higher  order  mode  groups,  because  higher 
order  mode  groups  have  a  larger  mode  count  and  carry  a  larger 
share  of  the  total  power  in  this  formulation.  In  modern  multi- 
mode  optical  communication  systems  where  laser  transmitter 
sources  are  commonly  used,  lower  mode  groups  often  carry  a 
larger  share  of  the  total  power  and  GIMFs  designed  based  on 
(14)  result  in  low  bandwidth. 

In  practice,  for  a  typical  untilted  center-launch  of  a  Gaussian- 
beam  laser,  zero-angular-momentum  modes  (m  =  0)  are  pre¬ 
dominantly  excited  and  mode  groups  g  =  1  and  g  =  3  carry 
the  largest  portion  of  the  power.  Therefore,  a  more  appropriate 
choice  of  aopt  can  be  be  expressed  as  =  P^l-^  or  even 
PgPg*  |g=i  =  0.  We  note  that  the  actual  value  of  a  in  high- 
bandwidth  commercial  GIMFs  is  proprietary  and  is  also  subject 
to  variations  in  the  manufacturing  process.  However,  because 
aopt  =  2  +  ?/  provides  a  better  modal  delay  balance  among  the 
lower  order  mode  groups  than  (14)  as  is  required  in  a  typical 
laser  launch,  we  speculate  that  it  is  a  more  accurate  representa¬ 
tion  of  the  refractive  index  profile  in  high-bandwidth  commer¬ 
cial  GIMFs  and  we  will  adopt  aopt  =2  +  2/  throughout  this 
paper. 


B.  Dispersion  Properties  of  a  Conventional  GIMF 

Conventional  high-bandwidth  commercial-grade  GIMFs 
such  as  Comings  InfiniCor  eSX+  are  commonly  optimized 
for  high-bandwidth  performance  at  850-nm  wavelength.  The 
typical  core  radius  R  =  25  pm  and  the  maximum  Ge-doping  is 
given  by  do  «  10.8.  At  A  =  850  nm,  this  results  in  a  =  2.103, 
no  =  1.4676,  A  =  0.01,  and  NA  =  0.208,  where  NA  is 
the  numerical  aperture.  We  reiterate  that  the  exact  values  of 
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the  aforementioned  parameters  are  often  proprietary,  but  it  is 
relatively  easy  to  speculate  on  the  values  of  these  parameters 
from  the  published  specifications  of  these  fibers  ( R  and  NA), 
and  the  preceding  analysis.  In  order  to  carry  out  specific  nu¬ 
merical  analysis  in  this  paper,  we  choose  the  aforementioned 
typical  parameters,  which  must  be  reasonably  close  to  those 
of  Comings  InfiniCor  eSX+  or  Thorlab’s  GIF50,  and  refer  to 
this  fiber  as  the  conventional  GIMF  (C-GIMF)  throughout  this 
paper.  All  the  numerical  results  in  this  paper  are  presented  for 
C-GIMF. 

In  Fig.  1,  we  plot  the  relative  delay  rg  \  =  fifjp  —  (delay 
relative  to  mode  group  g  =  1),  as  a  function  of  the  wavelength 
for  different  mode  groups.  We  observe  that  at  around  A  =  850 
nm  where  the  fiber  bandwidth  was  optimized  (by  choosing  a  = 
2+y  where  y  is  determined  at  A  =  850  nm),  the  relative  delay  is 
very  small  and  all  mode  groups  propagate  (approximately  )  with 
the  same  group  velocity.  However,  for  higher  mode  groups  such 
as  g  =  9.  t 9i  is  noticeably  larger  (0.51  fs/cm)  compared  with 
the  delay  experienced  by  lower  order  mode  groups  (e.g.,  T21  = 
0.018  fs/cm).  As  we  pointed  out  before,  such  high-order  mode 
groups  are  excited  less  efficiently  in  most  modern  multimode 
optical  communication  systems  and  optimal  values  of  a  should 
not  be  chosen  based  on  (14),  which  gives  a  high  weight  to  such 
high-order  mode  groups.  In  Fig.  1,  it  is  clear  that  at  wavelengths 
substantially  apart  from  A  =  850  nm,  the  modal  group  velocities 
are  not  longer  balanced,  but  are  “almost”  equally  spaced.  For 
example,  at  A  =  1.0  /jm,  T21  =  0.436  fs/cm  and  T31  =  0.896 
fs/cm. 

Earlier,  we  cautioned  that  the  presence  of  the  cladding  re¬ 
moves  the  degeneracy  of  the  delays  among  modes  within  each 
mode  group.  In  Fig.  2,  we  plot  the  relative  delay  at  A  =  850  nm 
for  all  the  modes  in  a  GIMF.  The  horizontal  axis  represents  the 
mode  group  number  g  and  the  dots  represent  the  values  of  rela¬ 
tive  delay  calculated  using  the  finite  element  method  (FEM)  of 
[17],  where  the  cladding  effect  is  now  taken  into  account.  The 
large  black  dots  are  for  zero-angular-momentum  modes  m  =  0, 
the  small  red  dots  belong  to  m  0  modes,  and  the  numerical 
values  of  the  delay  are  compared  with  the  analytical  expression 
from  (12)  plotted  in  the  green  dashed  line.  The  impact  of  the 
cladding  on  lifting  the  degeneracy  among  modes  in  higher  mode 
groups  is  more  pronounced,  as  expected.  In  Fig.  3,  we  make  a 
similar  plot,  except  at  A  =  1 000  nm,  where  we  observe  similar 
features.  As  expected,  the  delays  are  much  less  balanced  among 
the  mode  groups  for  A  =  1000  nm  compared  with  A  =  850  nm. 
In  this  paper,  as  we  will  discuss  later,  we  are  primarily  interested 
in  zero-angular-momentum  modes  (m  =  0),  and  Figs.  2  and  3 
clearly  illustrate  the  adequacy  of  analytical  expression  (12)  for 
g  <  11.  Figs.  2  and  3  are  quite  useful  in  estimating  the  nonlinear 
modal  interaction  length  between  different  mode  groups,  which 
is  limited  by  the  group  velocity  mismatch-induced  walk-off.  For 
example,  from  Fig.  3,  at  A  =  1000  nm,  the  relative  delay  be¬ 
tween  g  =  1  and  g  =  3  mode  groups  (consecutive  groups  with 
m  =  0  modes)  is  T31  =  0.896  fs/cm.  Therefore,  after  1  m  of 
propagation,  the  pulse  walk-off  is  around  90  fs,  resulting  in  very 
low  efficiency  interaction  between  pulses  with  duration  shorter 
than  100  fs.  However,  at  A  =  850  nm,  the  pulse  walk-off  is 
quite  smaller.  It  must  be  noted  that  these  are  crude  but  relevant 
approximations  and  one  must  include  such  effects  as  nonlinear 


trapping  in  order  to  obtain  a  more  reliable  estimate  of  the  true 
modal  nonlinear  interaction  length. 

Finally,  GVD  for  mode  group  g  =  1  (i.e.,  (3\  ')  is  plotted  as 
a  function  of  the  wavelength.  GVD  does  not  appreciably  differ 
among  mode  groups  and  is  taken  to  be  the  same  for  all  modes 
in  this  paper.  The  dashed  lines  identify  the  point  of  zero  GVD. 

C.  Approximate  Form  for  the  Propagation  Constant  (3 g 

It  is  sometimes  more  convenient  to  adopt  a  simplified  version 
of  (9),  when  possible.  For  A  <  1,  we  can  write 

Pg  ~  ri0k  ^1  -  Ax  -  ~Y~ )  +  0(A3).  (15) 


The  0(A2)  term  can  be  ignored  relative  to  the  O(A)  term  for 
short  propagation  lengths  where  it  does  not  amount  to  substan¬ 
tial  “extra”  phase  mismatch  between  different  mode  groups.  For 
example,  the  0(  A2)  phase  mismatch  between  g  =  1  and  g  =  3 
mode  groups  (consecutive  groups  with  m  =  0  modes)  is  given 
by 


603,1 


8A  L 
n0kR 2 


(16) 


where  a  «  2  is  assumed  for  simplicity  and  L  is  the  length  of 
the  GIMF.  Therefore,  the  0( A2)  term  contributes  a  2ir  addi¬ 
tional  phase  if  the  length  of  the  GIMF  is  L27r  =  nriokR2 /(AA). 
For  the  typical  design  parameters  of  conventional  commercial- 
grade  GIMFs  at  A  =  850  nm,  L27T  ~  54  cm;  therefore,  in  order 
to  reliably  ignore  the  0(  A2)  term,  the  length  of  the  GIMF  must 
be  considerably  shorter  than  L 2 -k-  We  note  that  is  even 
shorter  if  calculated  with  higher  order  mode  groups.  In  situa¬ 
tions  where  the  0(A2)  term  can  be  ignored  and  a  «  2  can  be 
trusted,  (15)  transforms  into  a  particularly  useful  form  of  linear 
in  g : 


R  ,  >/2A 

Pg  =  n0k0 - —g. 


111.  Pulse  Propagation  in  a  Nonlinear  MMF 

In  this  section,  we  briefly  review  the  basic  formulation  of  the 
nonlinear  pulse  propagation  in  a  MMF,  as  presented  by  Poletti 
and  Horak  [10],  Their  approach  is  based  on  the  pioneering  work 
of  Kolesik  and  Moloney  [18].  In  the  next  section,  we  will  adapt 
this  formulation  to  the  specific  case  of  the  nonlinear  pulse  prop¬ 
agation  in  a  GIMF. 

While  the  formulation  presented  in  [10]  includes  the  polar¬ 
ization  effects,  we  assume  in  this  paper  that  the  fields  are  ex¬ 
cited  and  remain  in  the  i-linear  polarization.  This  simplifying 
assumption  should  be  valid  based  on  our  experience  from  the 
nonlinear  pulse  propagation  in  SMFs  [5].  Essentially,  all  the 
nonpolarization-related  phenomena  in  fiber  optics  can  still  be 
observed.  While  the  polarization  effects  can  also  be  easily  in¬ 
cluded,  the  resulting  equation  is  unnecessarily  complicated  for 
the  purposes  of  this  paper,  which  is  focused  on  the  multi(spa- 
tial)-mode-related  behavior  of  of  the  nonlinear  propagation. 

The  modal  decomposition  of  the  (scalar)  electric  field  into  a 
sum  of  mode  function  Fp(x,y,uj)  with  envelopes  Ap(z,  cu)  is 
given  by 
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E(x,  y,z,t)  =  ^^2  j  dui  |  —ap(z,  u) 

x  exp  [i{f3p{(jj)z  —  wt)\  +  c.c}  (18) 

where  the  field  normalization  constant  Np(u>)  is  chosen  such 
that  the  mode  profile  satisfies  the  normalization  and  orthogo¬ 
nality  condition 

j  dxdyF*(x,y,uj)Fq(x,  y,w)  =  Spq  (19) 

and  \Ap(z,t)\2  represents  the  optical  power,  where  Ap(z,t)  is 
the  Fourier  transform  of  Ap(z,  u>). 

Poletti  and  Horak  [10]  have  shown  that  the  GNLSE  for  pulse 
propagation  in  multimode  fibers  can  be  written  as 

dAp(z,t) 

dz 


=  i  (AW  -  /3q0))  Ap{z,t) 


where  ujq  is  taken  as  the  carrier  frequency  of  the  initial  pulse. 
The  nonlinear  field  terms  are  given  by 

-P/mn  =  2 MZ,t)  J  dTR(T)Am{Z,t  ~  T)A*{Z,t-  t)  (21) 

Pi™n  =  J  dTR(r)Am{z,t  -  r)An{z,t  -  r)e2l“°T 

and  the  nonlinear  coupling  terms  are  given  by 

Q% L»  =  \J  dxdyF*(w)Fl(u)Fm(w)F*{u)  (22) 
QpilnJu)  =  ^  J  dxdyF*{w)F{,{w)Frn{w)Fn{w). 

The  nonlinear  response  function  can  be  written  as 

R(t)  =  (1  -  fR)8(r)  +  fRh(r )  (23) 


where  fR  «  0.18  is  the  fractional  contribution  of  the  Raman 
response  to  the  total  nonlinearity  and  h(r )  is  the  delayed  Raman 
response  function  [5],  The  shock  time  constants  are  given  by 


(b2)  ^ 

plmn 


(24) 


IV.  Nonlinear  Pulse  Propagation  in  a  GIMF 

The  GNLSE  (20)  is  our  starting  point  to  construct  a  suit¬ 
able  equation  for  the  nonlinear  pulse  propagation  in  a  GIMF.  A 


conventional  commercial-grade  GIMF  supports  more  than  100 
modes.  The  nonlinear  mode  propagation  can  turn  into  a  hopeless 
numerical  integration  of  more  than  100  coupled  NLSEs,  each 
with  more  than  one  million  nonlinear  coupling  terms.  While 
symmetry  arguments  can  somewhat  reduce  the  number  of  non¬ 
linear  coupling  terms  in  cylindrically  symmetric  optical  fibers  as 
shown  in  [10],  the  large  number  of  modes  and  coupling  terms 
can  easily  bury  the  interesting  physical  phenomena  and  espe¬ 
cially  make  it  difficult  to  use  the  observations  in  a  simple  form 
for  device  application;  therefore,  we  desire  to  restrict  our  study 
to  a  small  and  manageable  subspace  of  modes  in  a  consistent 
manner.  Fortunately,  not  only  is  this  possible  in  GIMFs,  it  is 
in  fact  a  practical  and  preferred  choice  [3].  In  this  paper,  we 
limit  our  formulation  and  studies  only  to  the  case  of  m  =  0  or 
zero-angular-momentum  modes.  This  choice  of  a  limited  sub¬ 
space  is  very  reasonable  and  desired,  because  the  injected  pulse 
is  often  coming  either  from  the  SMF,  which  is  fusion-spliced  to 
the  GIMF,  or  from  a  similar  azimuthally  symmetric  launch.  The 
experimental  evidence  presented  in  our  recent  work  in  [3]  con¬ 
firms  the  reliability  of  this  choice. 

It  is  very  important  to  note  that  in  the  nonlinear  case,  it  can  be 
shown  that  the  conservation  of  angular  momentum,  dictated  by 
the  symmetries  of  Qpl'mn  in  (20)  and  (22),  prevents  the  excita¬ 
tion  of  any  modes  with  m  ^  0  along  the  fiber  from  only  m  =  0 
modes;  therefore,  we  can  safely  limit  our  analysis  to  the  small 
and  manageable  zero-angular-momentum  subset  of  modes.  We 
note  that  m  ^  0  modes  may  still  be  excited  because  of  the  im¬ 
perfection  in  the  optical  fiber.  This  coupling  can  be  controlled 
and  minimized  in  short  fiber  segments.  Later  in  the  Appendix, 
we  will  briefly  discuss  some  of  the  issues  related  to  the  coupling 
among  modes  and  excitation  of  m  ^  0  modes,  especially  in  the 
context  of  macrobending. 

The  m  =  0  modes  of  the  GIMF  can  be  read  directly  from  (1), 
(3),  and  (4),  where  we  have 


Fp(uj)  =  Epfi(p,(/);uj) 


(25) 


Lp  is  the  Laguerre  polynomial,  which  is  simply  the  generalized 
Laguerre  polynomial  with  zero  angular  momentum.  After  a  few 
simple  algebraic  steps,  we  can  show  that 


Qplmn{w)  —  Qpllnni^)  —  Qp|m„(w) 


1  1 
3ttp5 


due  Lp  (u'jLj  (u)  Ljyi  (u)  . 

(26) 


From  (26),  we  can  readily  calculate  the  nonlinear  coupling 
among  m  —  0  modes  propagating  in  a  GIMF.  In  particular,  for 
the  lowest  order  mode,  we  obtain 

Qo00oH=3  2^f  =  3  4ff'  ^ 

We  note  that  pc  —  v/Spo  is  simply  the  modal  radius  of  this 
lowest  order  mode  and  is  the  key  spatial  dimension  in  a  GIMF 
that  sets  all  its  modal  properties,  not  the  core  radius  R  [3], 
Specifically,  pc  oc  \/A R  and  therefore,  the  effective  area  for 
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the  lowest  order  mode  A®s  =  7rp2  scales  as  the  core  radius,  not 
as  the  core  area  as  one  might  have  reasonably  imagined.  The 
scaling  of  the  modal  area  with  the  square  root  of  the  core  area 
has  important  consequences  for  nonlinear  pulse  propagation  in 
GIMFs.  For  example,  unlike  the  single-mode  fibers,  where  dou¬ 
bling  the  core  area  reduces  the  effective  nonlinearity  by  a  factor 
of  2,  the  effective  nonlinearity  of  each  mode  is  only  reduced  by 
a  factor  of  \/2  in  GIMFs. 

We  can  define  the  nonlinear  coupling  coefficient  in  GIMFs  as 

Iplmn  =  3  “  Qplmn(^)  ■  (28) 

We  note  that  the  forms  of  the  diagonal  terms  7 pppp  agree  with 
the  conventional  formulation  often  presented  for  the  SMF  and 
are  inversely  proportional  to  the  modal  area  of  that  particular 
mode  [5].  In  particular,  for  the  lowest  order  mode  (LP0i-like 
and  often  the  most  important  mode),  we  can  use  (26)-(28)  to 
show  that 

'«2^0 

7  =  7oooo  =  •  (29) 

^eft 


where  8(3 ^  =  (3^  -  (3 ^  and  8(3 [p^  =  (3^  -  (3^ .  A  final  and 
perhaps  major  simplification  is  possible  if  the  delayed  Raman 
response  and  the  shock  terms  are  ignored.  In  that  case,  the  last 
line  (nonlinear  coupling  terms)  in  (31)  can  be  simplified  as 

*7  E  PIplmnA-l^Z ,  t*)  An(z ,  f).  (32) 

/,m,n 

In  what  follows,  we  will  use  (3 1)  as  our  starting  point  to  analyze 
the  linear  and  nonlinear  propagation  of  optical  pulses  in  GIMFs. 
We  will  also  ignore  the  delayed  Raman  response  and  the  shock 
terms  and,  therefore,  adopt  the  approximation  presented  in  (32). 
The  impact  of  the  delayed  Raman  response  and  the  shock  terms 
will  be  addressed  elsewhere  in  the  future. 

V.  MMI  IN  GIMFs 

Perhaps,  the  most  trivial  application  of  the  GNLSE  (31)  is  to 
the  MMI  in  GIMFs  as  discussed  in  detail  in  [3]  and  [4].  We  con¬ 
sider  low-power  constant  wave  (CW)  light  propagation  where 
both  temporal  and  nonlinear  dynamics  can  be  ignored  and  ob¬ 
tain 


As  shown  in  (29),  we  will  simply  use  7  instead  of  70000  *n  the 
rest  of  this  paper. 

As  a  numerical  example,  we  consider  that  the  value  of  7  for 
a  typical  GIMF,  such  as  GIF50,  which  was  studied  in  [3],  is 
approximately  3.6  times  lower  than  that  of  a  typical  SMF-28, 
given  that  pc  «  7.73  /rm  for  GIF50  at  A  =  1550  nm.  The  value 
of  7  increases  at  shorter  wavelengths  and  becomes  larger  than 
1.2/W.km  at  A  «  1.0  ;/,m  (comparable  to  that  of  SMF-28  at 
A  =  1550  ran). 

Finally,  in  an  attempt  to  further  simplify  our  notation,  we  de¬ 
fined  rjpijnn  =  7 pimnh  with  70000  =  1-  We  note  that  r/p!mn  is 
a  fully  symmetric  tensor  and  represents  the  nonlinear  coupling 
among  propagating  modes  in  a  GIMF.  Here,  for  reference,  we 
calculate  the  nonlinear  coupling  terms  among  the  three  lowest 
order  modes  in  GIMFs 


_  i  _  _  3 

70001  —  7oon  —  7iin  —  2*  7ooi2  —  70022  —  g 

1  11 

70002  =  70111  —  70112  =  71112  =  Tjy  72222  =  7^ 


71122  =  — , 

lb 

11 

71222  - 


70122  -  ~ 
7,0222  =  h- 


(30) 


These  simplifications  allow  us  to  rewrite  the  GNLSE  (31)  in  a 
more  compact  form 

dAp(z,t) 

dz 

=  i8(3^Ap(z,t)-8(3^d-^^ 


*E  ^  (*flj  A-(M) 


(p) 


1>2 


d 


l,rri,n 


dt 

Vplrnn 


1  + 


ldt 


x  [p(+  +  P ^  1 1 
A  rlmn  '  rlm.n\  j 


(31) 


Ap(z)  —  Aj,(0)  exp  i8(3^z  .  (33) 


In  [3]  and  [4],  we  assumed  that  the  injected  beam  is  in  the  form 
of  a  Gaussian  with  radius  w,  centrally  aligned  with  the  GIMF, 
and  expressed  as 


H  = 


(34) 


where  (■ w\w )  =  1  and  the  bra-ket  notation  indicates  integration 
in  the  transverse  (p  —  <j>)  coordinates.  The  injected  beam  can  be 
easily  realized  from  the  nearly  Gaussian  mode  of  an  SMF  to  the 
GIMF.  The  initial  excitation  amplitudes  are  given  by  Ap(0)  = 
P'/'2  (Fp\w),  where  P0  is  the  power  of  the  input  Gaussian  field 
and  | Fp)  is  given  in  (25).  We  also  have 


{Fp\w) 


7+1  ’ 


$  = 


7-1 
7  +  1 ' 


(35) 


In  particular,  in  [3],  we  assumed  that  the  output  port  of  the  GIMF 
is  spliced  to  another  SMF  with  a  mode-field  radius  of  w' .  if  n 
is  the  power  carried  by  the  out-coupling  SMF,  the  power  cou¬ 
pling  efficiency  Pq/Pq  between  the  in-coupling  and  out-cou¬ 
pling  SMFs  can  be  written  as 


t(z) 


P' 

£_o 

Po 


p- 1 

J2W\Fp)(Fp\w)eis^z 

p= 0 


(36) 


Therefore,  we  get 


t(z)  = 


1 677' 


(7  +  l)2  (7/  +  l)2 


p-i 

y  -q/pqi'p 

p= 0 


eMp)* 


(37) 


For  m  =  0  modes  (g  =  2p  +  1),  from  (17)  we  have 

(3 p  =  n0k0  -  (2 p  +  l)ir/Z  8(3^  =  p/Z  (38) 


where  Z  =  nR/V2A.  If  we  use  (38)  in  (37),  we  obtain  the 
result  presented  in  (12)  of  [3]. 
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Fig.  4.  GVD  for  mode  group  g  =  1  is  plotted  as  a  function  of  the  wavelength. 
GVD  does  not  appreciably  differ  among  different  modes  and  can  usually  be 
taken  to  be  the  same  for  all  modes. 


Fig.  5 .  Relative  power  spectral  density  for  the  lowest  order  mode  in  GIMF  after 
0.1  m  of  propagation  in  the  C-GIMF.  The  injected  pulses  are  each  57  fs  long, 
with  200  kW  peak  power  at  Ap  =  805  nm  and  1  kW  peak  power  at  As  =  1030 
nm.  In  addition  to  the  expected  spectral  broadening,  the  FWM  generation  of  a 
new  peak  centered  at  A  a  =  660  nm  is  notable. 

VI.  New  Wavelength  Generation  Via  FWM 

In  the  previous  section,  we  used  the  GNLSE  (31)  to  explore 
a  MMI  problem  in  the  absence  of  both  dispersive  and  non¬ 
linear  effects  for  low-power  CW  light  propagation  in  GIMFs. 
In  this  section,  we  briefly  review  an  antipodal  application  of 
the  GNLSE  (31)  for  short  pulse  propagation  with  high  peak 
power.  In  particular,  we  will  observe  conversion  of  the  pump 
light  to  other  wavelengths  via  multimode  phase-matched  FWM 
in  GIMF. 

In  the  first  numerical  experiment,  we  launch  an  unchirped 
pulse  of  the  form  sech(f/T0)  with  T0  —  57  fs  and  peak  power 
Pq  —  200  kW  at  Xp  —  805  nm  in  a  C-GIMF  segment  in 
the  lowest  order  mode  (p.  m)  —  (0,  0);  this  pulse  can  be  ob¬ 
tained  from  a  typical  titanium-sapphire  mode-locked  laser.  A 
much  weaker  pulse  with  peak  power  1  kW  and  the  same  dura¬ 
tion  at  A s  —  1030  nm  is  also  launched  into  the  same  mode  of 
the  GIMF;  this  pulse  can  be  obtained  from  a  typical  Yb-doped 
mode-locked  fiber  laser.  Fig.  5  shows  the  “relative”  power  spec¬ 
tral  density  for  the  lowest  order  mode  after  propagating  a  length 
of  0.1  m  in  the  C-GIMF.  We  note  that  in  our  simulations,  we 
have  assumed  7  «  1.0/W.km  (typical  value  obtained  for  fused 
Silica  at  this  wavelength  from  (29)  for  C-GIMF),  which  results 
in  a  nonlinear  length  =  I/7P0  ~  5  mm.  In  addition  to  the 
expected  spectral  broadening  at  Ap  =  805  nm  due  to  dispersive 


Fig.  6.  Relative  power  spectral  density  for  the  lowest  order  mode  in  GIMF  after 
0.1  m  of  propagation  in  the  C-GIMF.  The  injected  pulses  are  each  5.6  ps  long, 
with  1  kW  peak  power  at  A s  =  1030  nm  and  1  W  peak  power  at  Ap  =  805 
nm.  The  FWM  generation  of  a  new  peak  centered  at  A^  =  660  nm  is  notable. 

and  nonlinear  SPM  effects,  the  generation  of  a  new  peak  cen¬ 
tered  at  A  ^  =  660  mu  is  notable.  The  conversion  of  light  from 
Ap  and  Xs  to  A^  is  mediated  by  a  FWM  process  where 

2  ujp  —  ws  +  wa  (39) 

where  ujp  =  2nc/Xp  and  so  on,  and 

2/3i(AP)  =  Pi(Xs)  +  P3(XA)  (40) 

where  (5\{Xp)  is  the  propagation  constant  of  the  lowest  order 
mode  at  the  pump  wavelength  Ap  (p  —  0  or  g  —  1  where 
according  to  (9)  g  is  used  as  a  lower  index  for  /3),  fii(Xs)  is 
the  propagation  constant  of  the  lowest  order  mode  at  the  Stokes 
wavelength  A s,  and  /33(A a)  is  the  propagation  constant  of  the 
next-order  mode  at  the  anti-Stokes  wavelength  A^  (p  =  1  gives 
g  =  3  for  the  index  of  /3).  We  note  that  the  FWM  process 
in  GIMFs  was  experimentally  observed  in  1981  by  Hill  et  al. 
[19];  using  the  formalism  presented  in  this  paper,  we  can  ex¬ 
plore  the  efficiency  of  such  processes  by  numerically  solving 
the  GNLSE  (31).  We  would  like  to  point  out  that  the  required 
phase  matching  presented  in  (40)  involves  different  modes  of 
the  GIMF,  as  also  observed  experimentally  in  [  1 9] .  A  more  com¬ 
prehensive  treatment  of  the  FWM  process  in  GIMFs  is  beyond 
the  scope  of  this  paper  and  will  be  presented  elsewhere.  The 
analysis  presented  in  this  section  is  merely  intended  to  show  the 
applicability  of  the  GNLSE  (31)  to  explore  the  diverse  and  rich 
linear,  dispersive,  and  nonlinear  phenomena  that  can  occur  due 
to  the  subtle  interaction  of  the  propagating  modes  in  GIMFs.  We 
note  that  all  the  features  observed  in  Fig.  5  (and  also  in  Fig.  6)  re¬ 
late  only  to  the  nonlinear  interactions  between  (m  =  0,p  =  0) 
and  (m  =  0,p  =  1)  modes. 

Before  we  close  this  section,  we  explore  a  similar  example, 
though  with  much  longer  pulses  and  lower  values  of  peak 
power.  Here,  we  launch  a  much  longer  unchirped  pulse  with 
T0  =  5.7  ps  and  peak  power  P0  —  1  kW  at  Ap  =  1030  nm 
and  a  much  weaker  pulse  with  peak  power  1  W  and  the  same 
duration  at  A s  —  805  nm,  both  into  the  lowest  order  mode  of 
the  GIMF.  Fig.  6  shows  the  “relative”  power  spectral  density 
for  the  lowest  order  mode  after  propagating  a  length  of  0. 1  m  in 
the  GIMF.  The  nonlinear  length  for  this  pulse  is  approximately 
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LNl  =  1  llPo  ~  1  m.  Considering  the  short  propagation  dis¬ 
tance  compared  with  the  nonlinear  length  and  the  long  duration 
of  the  pulse,  we  do  not  observe  any  appreciable  dispersive  or 
nonlinear  pulse  broadening  effects.  However,  the  impact  of 
the  FWM  processes  in  generating  new  wavelengths  can  be 
clearly  observed  in  Fig.  6.  We  note  that  other  phase-matched 
FWM  processes  beyond  (40)  are  possible  in  GIMFs  which  are 
responsible  for  the  generation  of  other  wavelengths  as  can  be 
also  observed  in  Fig.  6,  the  full  treatment  of  which  will  be 
presented  in  a  later  publication. 

VII.  Conclusion 

We  have  presented  a  detailed  account  of  the  modal  properties, 
dispersive  behavior,  and  nonlinear  mode  coupling  in  GIMF s  and 
derived  a  simplified  GNLSE,  which  can  be  used  to  explore  the 
rich  nonlinear  dynamics  related  to  the  propagation  and  inter¬ 
action  of  multiple  modes  in  GIMFs  in  a  tractable  manner.  Al¬ 
though  we  have  restricted  our  formulation  to  only  the  zero-an¬ 
gular-momentum  m  —  0  modes,  which  are  typically  less  than 
5  modes  among  the  more  than  100  modes  in  a  C-GIMF,  the  re¬ 
sulting  propagation  dynamics  is  still  considerably  more  com¬ 
plex  than  that  of  the  SMFs.  While  the  exclusion  of  m  ^  0  re¬ 
stricts  the  applicability  of  our  formalism  to  special  cases  and 
carefully  crafted  experiments,  we  anticipate  that  much  can  be 
learned  from  m  =  0  modes,  without  having  to  deal  with  a 
system  that  can  bury  much  of  the  interesting  physics  beneath 
its  overwhelming  complexity. 

We  note  that  this  paper  is  mainly  concerned  with  laying  out 
the  basic  formalism,  and  is  meant  to  pave  the  way  for  our  future 
studies  on  detailed  treatment  of  specific  nonlinear  phenomena 
such  as  SPM,  XPM,  and  FWM  in  GIMFs.  We  anticipate  that 
careful  and  detailed  studies  of  the  nonlinear  multimode  inter¬ 
actions  in  GIMFs  result  in  interesting  nonlinear  dynamics  and 
novel  device  applications. 


Appendix 

In  this  section,  we  will  present  order-of-magnitude  estimates 
of  mode  coupling  in  GIMFs;  in  particular,  we  will  discuss  the 
efficiency  of  m  ^  0  mode  excitations  from  m  =  0  modes. 

In  most  systems  involving  MMFs  and  especially  in  practical 
MMF  communication  systems,  linear  coupling  among  modes 
plays  an  important  role  and  must  be  properly  addressed.  The 
coupling  is  mediated  by  various  processes,  including  micro- 
and  macrobending,  stress,  cabling,  core  size  variations  along  the 
fiber,  and  various  other  manufacturing  defects.  Mode  coupling, 
because  of  its  random  nature,  is  traditionally  treated  stochasti¬ 
cally  in  MMF  communication  systems,  often  using  power-cou¬ 
pling  models  (see,  e.g.,  [20]— [23]).  As  a  starting  point  for  any 
method  to  address  the  mode  coupling  problem,  the  GNLSE  (20) 
must  be  supplemented  with  an  additional  mode  coupling  term 
i  CpqAq(z,  t )  on  the  right-hand  side,  where  Cpq  are  the  linear 
coupling  coefficients  among  the  modes,  representing  both  the 
deterministic  and  stochastic  processes.  It  must  be  noted  that  the 
presence  of  nonlinearity  often  complicates  the  stochastic  anal¬ 
ysis,  because  the  interaction  of  the  nonlinear  coupling  terms 
with  the  stochastic  linear  processes  may  result  in  colored  non- 
Gaussian  nonstationary  statistics  in  the  power  flow  equations. 
Unfortunately,  well-developed  methods  do  not  exist  to  handle 


such  stochastic  equations  and  one  must  resort  to  either  brute- 
force  application  of  Monte  Carlo  simulations,  or  exploring  the 
systems  in  the  low-power  regime,  where  the  nonlinear  terms  can 
be  linearized. 

For  the  purposes  of  this  paper,  which  is  more  concerned  with 
short  segments  of  GIMFs,  macrobending  is  the  most  important 
coupling  process.  In  particular,  macrobending  can  potentially 
result  in  the  excitation  of  nonzero-angular-momentum  modes 
(to  =  ±1)  from  m  =  0  modes.  The  coupling  between  (to  = 
0 ,p  =  0)  and  (to  =  \,p  =  0)  modes  can  be  easily  estimated 
using  the  formalism  presented  in  [1]  and  [24],  and  can  be  written 
as 


C 


npkp 

Rr 


(Ira:  0p|ir|0m!  Op) 


noirpo 

Rc\ 


(41) 


where  Rc  is  bending  radius  of  the  GIMF.  |0m,  0p)  and  |lm,  0p) 
represent  the  (to  =  Q,p  =  0)  and  (to  =  1  ,p  =  0)  modes,  re¬ 
spectively.  The  length  scale  associated  with  the  coupling  term 
is  Lc  =  1/C;  for  the  design  parameters  of  C-GIMF  at  A  =  1550 
nm,  Cc  «  RJ 16.  In  practice,  the  coupling  between  these  two 
modes  is  strongly  suppressed  due  to  the  large  phase  mismatch, 
because  the  two  modes  belong  to  two  different  mode  groups. 
The  phase  mismatch  A/3  between  these  modes  can  be  approxi¬ 
mated  using  (17),  resulting  in  A/3  =  \/2A/R.  The  power  cou¬ 
pling  efficiency  (the  maximum  fraction  of  the  power  in  to  =  0 
mode  that  couples  to  the  to  =  1  mode)  can  be  characterized  by 
(  =  C2/(  A/3/2)2,  where  it  can  be  easily  shown  that  (-Cl, 
unless  Rc  <  6  mm. 

From  the  aforementioned  analysis,  it  is  clear  that  for  device 
applications  which  rely  on  short  segments  of  GIMFs,  and  when 
only  m  =  0  modes  are  initially  excited  as  in  [3],  [4],  excitation 
of  to  =  ±1  modes  is  a  highly  inefficient  process;  for  example, 
for  a  fiber  segment  with  R,c  « 1  m,  (  «  3  x  10~5 .  As  long  as  the 
fiber  is  held  fairly  straight,  the  power  coupled  to  the  to  =  ±1 
mode  due  to  macrobending  is  negligible.  We  emphasize  that 
such  a  large  phase-mismatch  between  to  =  0  and  to  =  ±1 
modes  is  a  characteristic  of  GIMFs  and  does  not  hold  in  general 
for  all  MMFs;  this  is  one  of  the  key  reasons  in  making  GIMFs 
so  attractive  for  device  applications  [3],  [4]. 

Core  size  distortions  can  also  lead  to  coupling  among  modes 
with  Am  =  ±2  [25],  Such  distortions  are  the  result  of  manufac¬ 
turing  errors  and  are  often  small;  however,  A  to  =  ±2  processes 
can  result  in  couplings  within  a  mode  group  where  the  modes  are 
approximately  phase-matched.  Detailed  analysis  of  this  issue  is 
beyond  the  scope  of  this  paper  and  the  results  depend  strongly 
on  the  manufacturing  quality  of  the  fiber.  For  short-length  de¬ 
vice  applications  such  as  those  in  [3]  and  [4],  this  effect  is  neg¬ 
ligible,  for  the  range  of  the  fibers  studied  in  [3]  and  [4].  How¬ 
ever,  in  MMF  communication  systems  where  propagation  dis¬ 
tances  are  much  longer,  one  must  worry  about  both  A  to  =  ±1 
and  A  to  =  ±2  process,  where  partial  phase  matching  in  the 
A  to  =  ±1  processes  can  also  be  assisted  by  various  random 
variations  along  the  fiber. 
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1.  INTRODUCTION 

Multimode  interference  (MMI)  in  optical  fibers  has  been  used 
successfully  in  recent  years  for  various  device  applications, 
including  beam  shapers,  sensors,  and  filters  [1-8].  MMI  in  a 
graded-index  multimode  optical  fiber  (GIMF)  was  recently 
used  to  create  very  low-loss  couplers  between  two  single¬ 
mode  optical  fibers  (SMFs)  with  very  different  mode-field 
diameters  [9,H)].  Here,  the  linear  analysis  of  [9]  is  extended 
and  the  MMI  effect  in  the  nonlinear  regime  for  the  GIMF  cou¬ 
pler  geometry  shown  in  Fig.  J_  is  investigated.  The  main  inten¬ 
tion  of  this  study  is  to  explore  the  possibility  of  using  the 
simple  SMF-GIMF-SMF  geometry  in  nonlinear  device  applica¬ 
tions,  such  as  in  optical  signal  processing  [11]  or  as  a  saturable 
absorber  (SA)  in  mode-locked  fiber  lasers  [12, 13], 

GIMFs  are  commonly  used  in  fiber-optic  communications 
to  reduce  modal  dispersion,  as  the  group  velocities  of  all 
modes  are  nearly  identical  at  the  design  wavelength  [14,15], 
However,  GIMFs  exhibit  another  unique  property  that  makes 
them  very  attractive  for  MMI  applications:  the  propagation 
constants  of  their  modes  are  equally  spaced.  Consequently, 
their  self-imaging  lengths  can  be  very  short,  even  less  than 
1  mm  [9];  therefore,  it  is  possible  to  make  extremely  short 
(even  submillimeter)  practical  MMI  devices  [10]  or  to  easily 
tune  the  length  of  the  GIMF  coupler  for  a  specific  MMI 
application  [10], 

In  this  paper,  our  studies  are  focused  on  the  simplest  case 
of  nonlinear  MMI  (NL-MMI)  in  the  setup  shown  in  Fig.  1_,  with 
identical  input  and  output  SMFs.  This  choice  is  also  the  most 
practical  one  for  most  device  applications.  Our  results  will 
elucidate  some  general  behaviors  of  such  couplers  in  the  non¬ 
linear  regime.  More  complex  geometries  involving  different 


fiber  junctions  and  using  various  mode  conversion  techniques 
[16],  such  as  long-period  gratings  [17],  can  modify  some  of  the 
observations  and  conclusions.  However,  such  complex  mod¬ 
ifications  are  less  likely  to  be  adopted  in  practical  device  ap¬ 
plications  in  the  near  future.  We  hope  that  our  results  can  be  a 
useful  starting  point  for  future  studies  of  NL-MMI  and  mode 
conversion  in  more  complex  systems. 

NL-MMI  has  been  studied  in  various  contexts  over  the  years 
(e.g.,  [18]).  Generally  speaking,  nonlinear  devices  that  operate 
based  on  nonlinear  mode  switching  and  coupling  in  spatially 
separated  waveguides  can  also  be  viewed  as  MMI  devices 
[19-22],  Nonlinear  polarization  rotation  is  another  important 
example  of  NL-MMI  (between  two  orthogonal  polarization 
modes)  that  closely  resembles  our  analysis  here  [23,24], 

In  Section  2,  a  general  overview  of  GIMFs,  nonlinear  propa¬ 
gation  in  multimode  optical  fibers,  and  MMI  are  presented. 
NL-MMI  in  GIMFs  in  the  context  of  the  SMF-GIMF-SMF  geom¬ 
etry  of  Fig.  1  is  discussed  in  Section  3.  In  Subsection  3.A,  a 
reduced  version  of  the  model  that  includes  only  two  propagat¬ 
ing  modes  is  analyzed.  This  simplification  helps  us  to  develop 
proper  insight  into  the  physics  of  NL-MMI  in  GIMFs.  In  Sub¬ 
section  3.B,  our  analysis  is  extended  to  the  realistic  case  of 
five  propagating  modes.  The  impact  of  the  presence  of  addi¬ 
tional  propagating  modes  in  the  NL-MMI  setup  can  be  shown 
by  comparing  these  results  with  the  results  of  Subsection  3.A. 
In  Section  4,  the  saturable  absorption  behavior  of  the  SMF- 
GIMF-SMF  geometry  is  analyzed  as  an  example  of  the  utility 
of  this  setup  for  device  applications.  The  concluding  remarks 
are  presented  in  Section  5,  where  it  is  shown  that  the  SMF- 
GIMF-SMF  geometry  can  be  a  viable  design  for  nonlinear 
switching  or  saturable  absorption. 
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Fig.  1.  GIMF  of  length  L  is  used  as  an  intermediate  coupler  between 
two  SMF  fibers.  In  [9,10],  this  geometry  was  used  to  create  very 
low-loss  couplers  between  two  SMFs  with  very  different  mode-field 
diameters.  In  this  paper,  nonlinear  MMI  effects  for  identical  input 
and  output  SMFs  are  explored. 
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The  formalism  and  results  in  this  paper  are  mainly  laid  out 
in  dimensionless  units.  This  choice  is  common  in  nonlinear 
fiber  optics  because  it  reduces  the  number  of  parameters 
and  prevents  redundancy  in  the  analysis.  For  real-world  appli¬ 
cations,  it  is  easy  to  convert  back  to  the  dimensionful  param¬ 
eters.  In  Section  5,  the  performance  of  a  NL-MMI  SA  device 
with  specific  parameters  from  commercially  available  optical 
fibers  is  analyzed. 

In  order  to  reduce  the  complexity  of  the  analysis  and  to 
make  the  problem  more  tractable,  the  studies  in  this  paper 
are  considered  in  the  continuous  wave  (CW)  limit.  In  practice, 
after  the  design  parameters  are  chosen  based  on  the  CW 
analysis,  the  temporal  effects  can  be  included  to  optimize 
the  final  design.  GIMFs,  if  used  near  their  optimal  design 
wavelengths,  are  particularly  attractive  for  short-pulse  ultra¬ 
fast  applications  compared  with  other  MMFs.  The  low  modal 
dispersion  in  GIMFs  ensures  that  pulses  do  not  break  up  in 
short  GIMF  segments  [26],  making  the  CW  analysis  an 
adequate  approximation. 

Finally,  it  should  be  pointed  out  that  the  transmission 
through  the  SMF-GIMF-SMF  geometry  is  a  periodic  fimction 
of  the  frequency  of  the  light  source  [9,10].  For  mode-locking 
applications  where  short  pulses  are  generated  with  very  large 
spectral  bandwidth,  the  length  of  the  GIMF  needs  to  be  suffi¬ 
ciently  short  to  provide  the  necessary  spectral  transmission 
window.  Moreover,  the  SMF-GIMF-SMF  geometry  can  be  po¬ 
tentially  used  as  the  spectral  filter  that  is  required  to  stabilize 
the  mode-locking  operation  of  normal-dispersion  high-energy 
femtosecond  fiber  lasers  [26] . 

2.  FUNDAMENTALS 

In  the  following  three  subsections,  a  brief  general  overview  of 
GIMFs  (Subsection  2. A),  nonlinear  multimode  propagation  of 
light  in  GIMFs  (Subsection  2.B),  and  the  general  formulation 
of  the  MMI  phoenomenon  in  GIMFs  in  the  context  of  the  SMF- 
GIMF-SMF  coupler  (Subsection  2.C)  will  be  presented.  The 
formulation  presented  in  this  section  will  be  used  to  analyze 
the  nonlinear  behavior  of  the  SMF-GIMF-SMF  geometry  in 
detail  in  subsequent  sections. 


EP.m(P-  <t>)=N%-^  exp  (-  I^Ljr'  (2) 

where  L\T  are  generalized  Laguerre  polynomials,  and  p0  and 
N™  are  given  by 


Po  ~  (fc0n0)1/2( 2A)1/4  ’  P  _  V  *(p  +  M)f  1  J 

where  k0  =  2n/X.  The  coefficients  Np  of  these  Laguerre- 
Gauss  modes  (LG pm)  are  chosen  such  that  the  modes  are  or¬ 
thonormal.  Using  the  bra-ket  notation,  the  electric  field  profile 
of  the  LGpm  mode  is  identified  as  \p.  to)  =  Epm(p ,  <p)  and  the 
orthonormality  condition  is  expressed  as  ( p,m\p',m ')  = 
Sp  p'Sm,m',  where  the  bra-ket  indicates  integration  in  the  trans¬ 
verse  (f >-(/>)  coordinates. 

All  modes  with  equal  group  mode  number  g  =  2 p  +  \m\  +  1 
are  almost  degenerate  in  the  value  of  the  propagation  con¬ 
stant,  which  is  given  according  to  the  formula 

fiun  =  n0k(l-2AXg)1'2,  (4) 


where 


/  g  \2a/(a+2) 

= 


(5) 


Na  is  the  total  number  of  guided  modes  (counting  the  polari¬ 
zation  degeneracy)  and  is  given  by 


Na  =  ^n20klR*A. 


(6) 


While  the  analysis  presented  in  this  article  is  generally 
applicable  to  a  wide  range  of  GIMFs  and  is  presented  in  di¬ 
mensionless  units,  for  specific  numerical  arguments  and  com¬ 
parisons,  the  parameters  of  a  conventional  high-bandwidth 
commercial-grade  GIMF  such  as  Coming’s  InfiniCor  eSX+  will 
be  used;  this  fiber  is  optimized  for  high  bandwidth  perfor¬ 
mance  at  850  nm  wavelength,  with  the  core  radius  of 
R  =  25  pm.  The  specified  GIMF  will  be  referred  to  as  C-GIMF 
standing  for  the  conventional  GIMF.  It  should  be  noted  that 
the  key  spatial  dimension  in  a  GIMF  that  sets  all  its  modal 
properties  is  not  the  core  radius  R,  but  is  pc  =  V2p0,  where 
Po  is  defined  in  Eq.  (3).  pc  is  the  mode  radius  of  the  LG0o  mode; 
for  C-GIMF,  pc  a  7.7  |im  at  1550  nm  wavelength  and  pc  « 
5.7  (im  at  850  nm  wavelength. 


A.  Overview  of  GIMFs 

The  refractive  index  profile  of  a  GIMF  is  given  by 


n2(p)  =  Wq 


(1) 


where  R  is  the  core  radius,  n0  is  the  refractive  index  at  the 
center  of  the  core,  A  is  the  index  step,  a  «  2  characterizes  a 
near  parabolic-index  profile  in  the  core  (p  <  R),  and  a  =  0 
in  the  cladding  (p  >  R).  The  transverse  electric  field  profile 
of  a  confined  mode,  with  radial  p  (p  >  0)  and  angular  to 
integer  numbers,  can  be  expressed  as  [25] 


B.  Overview  of  Nonlinear  Propagation  in  GIMFs 

Consider  an  input  optical  field  injected  into  the  core  of  the 
GIMF.  The  injected  field  excites  the  LGpm  modes  of  the  GIMF 
with  different  amplitudes  A;,m(0).  Apm(z)  will  be  regarded  as 
the  envelope  of  the  electric  field  and  0  =  0  indicates  the  input 
end  of  the  GIMF  [25],  In  order  to  simplify  our  notation,  the 
collective  index  p  will  be  used  to  represent  the  p,  m  index  pair 
of  the  LGpm  modes,  i.e.,  g  =  (p.  m)  in  the  rest  of  this  section. 
However,  when  necessary,  we  will  directly  use  the  (p,  to)  la¬ 
bels  or  specify  how  to  translate  the  p  index  into  the  (p,  to) 
index  pair.  \A„(z)\2  represents  the  optical  power  in  the 
LGpm  mode. 
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In  the  scalar  approximation,  the  generalized  nonlinear 
Schrodinger  equation  (GNLSE)  describing  the  CW  longi¬ 
tudinal  evolution  of  Apm(z)  can  be  written  as 

^  =  i5PpAp  +  v/Yjl^K^AAt  ■  (7) 

V,K,£ 

where  5pp  =  Pp  -  pp  is  the  propagation  constant  of  the 
mode  with  the  collective  index  //  and  /?(1)  is  the  propagation 
constant  of  the  LG00  mode  with  g  =  1  from  Eq.  (4).  The  “nor¬ 
malized”  nonlinear  coupling  coefficient  is  a  fully  symmetric 
tensor  and  is  defined  as  tjpl/K^  =  y^/y  where 

W  =  (n;J°)  /  d2xE*EvEKEf.  (8) 


defined  as  |in),  where  (in|in)  =  1,  the  mode  amplitude  in  the 
GIMF  at  the  input  can  be  written  as  A  ( 0)  = 
where  \p.  to)  is  the  GIMF  mode  defined  in  Eq.  (2). 

The  relative  transmitted  power  to  the  output  SMF  is 
given  by 

T  =  l\(out\E(p,<f>,L))\\  (13) 

where  |out)  is  the  normalized  mode  of  the  output  SMF,  and 
| E(p,<p,L))  is  the  total  field  at  the  output  facet  of  the  GIMF 
with  a  length  of  L,  calculated  from 

E(p,  (p,L)  =  eW^M^E/p,  #).  (14) 


In  Eq.  (8),  Ep  is  the  shorthand  notation  for  Epm(p,  if>)-y  =  yoooo 
is  the  nonlinear  coefficient  of  the  LG00  mode,  n2  is  the  non¬ 
linear  index  coefficient,  and  a>0  is  the  carrier  frequency. 
For  the  lowest  order  mode  LG0o  for  which  pc  is  simply  the 
modal  radius,  one  can  define  y0ooo  =  fcffloMjj,  where 
=  np%c.  The  total  optical  power  in  the  GIMF  is  given  by 
P  =  and  is  conserved  in  propagation  along  the 

GIMF,  i.e.,  d,P  =  0. 

It  is  more  convenient  to  express  the  GNLSE  in  dimension¬ 
less  units.  /:(,  can  be  defined  as 

Bf,(z)  =  (9) 

VP 

and  rescale  the  longitudinal  coordinate  z  by  the  difference  be¬ 
tween  the  propagation  constants  of  the  first  and  the  second 
mode  groups  pm  -  /?(2), 


Identical  input  and  output  SMFs  are  the  focus  of  this  paper 
(|in)  =  |out)),  for  which  the  relative  power  transmission  can 
be  written  as 


T 


1 

p2 


£/*(0)A„(*) 


!>*('  0)Bp(z) 


(15) 


The  injected  beam  from  the  input  SMF  is  in  the  form  of  a 
Gaussian  with  radius  w  and  can  be  expressed  as 

I  in)  =  exp  (-  P  \ .  (16) 

V  nwi  \  wr ) 

Because  of  the  azimuthal  symmetry  of  the  input  Gaussian 
beam,  if  the  input  SMF  is  centrally  aligned  with  the  GIMF,  only 
LGpo  (to,  =  0)  modes  can  be  excited  in  the  GIMF.  The  excita¬ 
tion  amplitude  can  be  calculated  as  [9] 


C  =  zx(pm-p(2)).  (10) 

Using  these  transformations,  Eq.  (7)  can  be  simplified  as 


Bp(0)  = 


2>  q/p 
>1+1 


y  =  V-1 
>1  +  1 ' 


(17) 


=  -i(rH  +  y)Bp  +  if  ^  ,  (11) 

where  Eq.  (11)  is  expressed  in  terms  of  the  dimensionless 
coefficients 

r  =  P?-  Pm 

h  P(  2)  -  Pm' 

and  can  be  solved  for  the  dimensionless  field  Bfl  in  tenns  of 
the  dimensionless  longitudinal  coordinate  t .  It  can  be  realized 
that  the  rescaled  fields  Bfl  satisfy  the  power  conservation  con¬ 
dition  J2fi\Bf,(z)\2  =  1  at  any  points  along  the  GIMF. 

Using  Eqs.  (4)  and  (5)  the  fact  that  Xg  -c  1,  one  can  show 
that  rfl  is  non-negative  and  is  nearly  an  integer:  rfl  agp-  1, 
where  gfl  is  the  group  number  associated  with  mode  p  [10], 

C.  Overview  of  MMI  in  GIMFs 

The  particular  setup  that  is  considered  in  this  paper  is  shown 
in  Fig.  1  and  consists  of  injecting  a  nearly  Gaussian  beam  from 
the  input  SMF,  which  is  spliced  to  the  input  facet  of  the  GIMF, 
and  collecting  the  light  at  the  other  end  of  the  GIMF  from  the 
output  SMF.  If  the  normalized  mode  of  the  input  SMF  is 


In  Fig.  2,  the  excitation  amplitudes  of  the  LGp0  modes  are 
displayed  as  a  function  of  the  radial  number  p,  for  three  differ¬ 
ent  values  of  the  ij  parameter.  For  the  case  of  p  =  1  (not 
shown  in  Fig.  2),  only  the  LG0o  mode  is  excited,  because 
the  input  Gaussian  beam  is  mode-matched  to  the  LG0o  mode 
of  the  GIMF.  For  a  slightly  larger  value  of ;/  =  1.1,  only  0.2%  of 
power  is  coupled  to  the  LG0i  mode  and  the  power  coupled  to 
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Fig.  2.  Excitation  amplitudes  of  the  LGp0  modes  from  Eq.  (17)  are 
plotted  as  a  function  of  the  radial  number  p,  for  three  different  values 
of  the  parameter. 


Y  = 


yP 


Pm  ~  Pm 


(12) 
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higher  order  modes  is  negligible.  For  the  larger  value  of ;/  =  4, 
only  64%  of  the  power  is  coupled  to  the  LG0o  mode.  In  this 
case,  the  relative  power  coupled  to  the  LG01,  LG02,  and 
LG03  modes  are  23%,  8%,  and  3%,  respectively.  //  =  4  is  ob¬ 
tained  when  coupling  the  C-GIMF  to  an  SMF  with  the 
mode-field  diameter  of  nearly  7.7  pm  at  1560  nrn  wavelength. 
For  coupling  of  the  C-GIMF  to  Coming  SMF-28TM  (SMF-28) 
with  the  mode-field  diameter  of  10.4  pm  at  1550  nm  wave¬ 
length,  the  value  of  ij  is  nearly  equal  to  2.2. 

There  are  two  reasons  for  not  exploring  the  case  of  ij  <  1  in 
this  manuscript.  First,  under  the  rj  -*  r) -I  transformation, 
Eq.  (17)  remains  unchanged  except  for  an  overall  multiplica¬ 
tive  factor  of  (-l)3*  [9],  Despite  this  change  in  the  initial  phase, 
no  new  insight  is  obtained  in  overall  observations  and  general 
conclusions  presented  in  this  paper.  Second,  ;/  is  larger  than 
unity  for  most  combinations  of  commercially  available  SMFs 
and  GIMFs.  Therefore,  our  analysis  is  limited  to  ij  >  1  to  avoid 
adding  unnecessary  complexity. 

3.  NL-MMI:  SMF-GIMF -SMF  GEOMETRY 

In  the  following,  our  studies  will  only  include  the  case  of  m  = 
0  or  zero-angular  momentum  modes  (see  the  discussion  in 
Subsection  2.C).  The  experimental  evidence  presented  in 
[9,10]  confirms  the  reliability  of  this  choice.  Moreover,  using 
the  formalism  outlined  in  [25]  it  can  be  shown  that  the  con¬ 
servation  of  angular  momentum,  dictated  by  the  symmetries 
of  the  nonlinear  coupling  terms,  prevents  the  excitation  of  any 
modes  with  m  ^  0  along  the  fiber  from  only  m  =  0  modes. 
Therefore,  our  analysis  can  be  safely  limited  to  the  small 
and  manageable  zero  angular  momentum  subset  of  modes. 
It  is  worth  mentioning  that  m  ^  0  modes  may  still  be  excited 
because  of  the  bending  and  imperfections  in  the  fiber.  This 
issue  will  be  addressed  in  Section  5,  but  for  now  our  analysis 
is  limited  to  the  subspace  of  m  =  0  modes. 

Because  only  LG?)0  modes  are  considered  in  the  analysis, 
one  can  readily  identify  the  collective  index  p  used  in  Section  2 
with  the  radial  mode  number  p. 

Therefore,  all  the  indexes  in  the  relevant  equations  such  as 
in  Eqs.  (7),  (8),  (11),  and  (14)  can  be  directly  replaced  with  the 
radial  mode  numbers  of  the  LGf)0  modes. 

In  a  typical  GIMF,  only  a  handful  of  m  =  0  modes  are  sup¬ 
ported  in  the  core;  for  example,  only  five  m  =  0  modes  are 
supported  at  1550  nm  wavelength  in  a  C-GIMF.  Even  for  such 
a  small  number  of  modes,  the  GNLSE  will  have  many  terms 
and  becomes  quite  complicated.  In  order  to  develop  the 
proper  fundamental  understanding  and  intuition  on  the  key 
mechanisms  involved  in  the  NL-MMI  of  these  modes  in  a 
GIMF,  our  analysis  is  initially  limited  to  the  smaller  subset 
of  LG0o  and  LG10  modes  in  Subsection  3. A.  This  way,  the  phys¬ 
ics  of  the  NL-MMI  in  the  SMF-GIMF-SMF  is  not  buried  under 
the  complexity  introduced  by  the  large  number  of  modes 
propagating  in  the  GIMF.  Eventually,  in  Subsection  3.B,  the 
number  of  propagating  modes  in  the  GIMF  will  be  increased 
to  five,  considering  LG00  through  LG40,  and  the  similarities  and 
differences  between  the  two-mode  and  the  five-mode  scenar¬ 
ios  will  be  explored. 

As  a  reference  for  comparison  with  nonlinear  couplings  for 
self-phase  modulation  and  cross-phase  modulation  in  conven¬ 
tional  SMFs  [24],  the  nonlinear  coupling  terms  among  the 
three  lowest  order  modes  in  GIMFs  (i.e.,  LG00,  LG10,  and 
LG2o),  are:  7/0ooi  =  *7oon  =  *7im  =  1/2,  */ooi2  =  t/0022  =  3/8, 


*72222  —  H/32,  1)1122  —  5/16,  7/0002  —  '7oin  —  *7oii2  —  *71112  — 

1/4,  7/0122  =  3/16,  7/1222  =  H/64,  7/0222  =  5/32. 

A.  NL-MMI  with  and  Modes 

In  this  subsection,  it  is  assumed  that  only  LG0o  and  LG10 
modes  are  excited  in  the  GIMF  and  the  optical  power  is  dis¬ 
tributed  only  between  these  two  modes.  It  should  be  taken 
into  account  that  this  assumption  is  only  strictly  valid  for  1/  « 
1  from  Eq.  (17)  (p0  »fh  in  that  case).  However,  our  intention 
in  this  section  is  not  to  simulate  an  exact  experimentally  real¬ 
izable  configuration  in  the  spirit  of  Fig.  1.  Rather,  this  limited 
two-mode  subspace  is  used  to  gain  insight  into  the  physics  of 
NL-MMI  in  GIMFs.  Therefore,  in  this  subsection,  the  limita¬ 
tions  of  Eq.  (17)  are  abandoned  and  it  is  assumed  that 
p0  +  Pi  =  1,  where  p0  =  |fi0|2  and  Pi  =  I-B1I2  can  adopt  any 
(positive)  values  subject  to  this  condition,  hi  practice,  it 
may  be  possible  to  create  this  scenario  experimentally  in  a 
more  complex  setup  than  that  of  Fig.  1,  e.g.,  by  using  mode 
conversion  techniques.  However,  such  details  are  beyond 
the  scope  and  intentions  of  this  subsection. 

Using  Eq.  (11),  the  following  coupled  nonlinear  equations 
are  obtained: 

d^B0  =  iy(\B^  +  \\Bl\^Bl  +  |  (B20B*  +  B'fB* ) , 

d(B,  =  -inBj,  +  f  A)  +  f  (-l-Bil2#!  +  BIB*  +  \B\B*y 

(18) 

Nonlinear  effects  can  be  ignored  when  /  L  «0,  where  L  =  Lx 
(flnj  -  /7(2j)  and  L  is  the  total  length  of  the  GIMF.  In  the 
absence  of  nonlinear  effects,  the  solution  to  Eq.  (18)  can 
be  written  as 

B0(O  =  fib(0),  Bi(0  =  e-<r‘ABi(0).  (19) 

In  the  linear  case,  the  relative  power  of  each  mode  is  con¬ 
served,  i.e.,  p0  and  p\  remain  unchanged  along  the  GIMF. 
The  relative  power  transmission  in  an  SMF-GIMF-SMF  setup 
in  the  linear  limit  with  two  modes  can  be  calculated  from 
Eq.  (15)  as 

t  =  1  -  4p0Pi  sin2  ^0 .  (20) 

The  relative  power  transmission  is  a  periodic  function  of  the 
length  of  the  GIMF  and  varies  periodically  between  rIllin  =  1  - 
4p0Pi  and  rmax  =  1;  rmin  =  0  in  the  special  case  where  the 
power  is  equally  distributed  (p0  =  pi  =  1/2).  These  issues 
have  been  discussed  in  further  detail  in  [9],  where  the  linear 
MMI  in  GIMFs  is  explored  in  the  presence  of  multiple 
m  =  0  modes. 

The  main  parameters  that  determine  the  value  of  the  rela¬ 
tive  power  transmission  in  the  SMF-GIMF-SMF  setup  are  the 
relative  power  of  the  modes  (p0  and  pt  ),  the  normalized  GIMF 
length  (L),  and  the  normalized  nonlinear  coefficient  y.  In  Fig.  3, 
the  relative  power  transmission  is  plotted  as  a  function  of  L 
for  three  different  cases  of)  =  0  (solid),  y  =  0.7  (dashed),  and 
y  =  3  (dotted).  All  cases  are  plotted  for  relative  power  values 
of  p0  =  pi  =  0.5.  The  linear  case  with  y  =  0  clearly  follows 
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solid,  dashed,  dotted:  f= 0,  0.7,  3 


Fig.  3.  Relative  power  transmission  is  plotted  as  a  function  of  the 
normalized  GIMF  length  L  for  the  case  of  LG00  and  LG10  modes  when 
Pa  =  Pi  =  0-5  at  f  =  0,  for  y  =  0  (solid),  y  =  0.7  (dashed),  and  y  =  3 
(dotted). 

Eq.  (20),  where  the  oscillation  periodicity  of  r  with  respect  to 
L  (normalized  GIMF  length)  is  T L  =  k.  k. 

When  the  nonlinear  effects  are  small,  the  nonlinearity 
merely  results  in  an  additional  phase  mismatch  between 
the  LG00  and  LG10  modes,  e.g.,  for  y  =  0.7  the  accumulated 
nonlinear  phase  mismatch  is  nearly  n  at  L  =  4.5/r  and  is 
responsible  for  changing  the  r  =  0  in  the  linear  case  to 
1  for  y  =  0.7  in  Fig.  3.  It  should  be  emphasized  that  in  the  pres¬ 
ence  of  nonlinearity,  p0  and  p1  are  not  conserved  and  vary 
periodically  along  the  GIMF  (LG0o  and  LG10  exchange  power). 
However,  unless  we  explicitly  specify  their  f -dependence, 
when  referring  to  p0  and  p,,  we  mean  their  initial  value  at 
the  f  =  0  point  in  the  GIMF  and  refrain  from  introducing 
new  variables. 

The  linear  and  nonlinear  cases  of  y  =  0  and  y  =  0.7  differ  in 
another  important  feature  besides  the  cumulative  nonlinear 
phase  mismatch;  the  minimum  value  of  r  for  y  =  0.7  is  larger 
than  zero.  This  effect  becomes  more  prominent  for  larger  rel¬ 
ative  nonlinear  coefficients.  In  fact,  as  can  be  seen  for  y  =  3  in 
Fig.  3,  the  relative  power  transmission  remains  above  60%  for 
all  values  of  L. 

The  relative  power  transmission  explored  in  Fig.  3  is  for  the 
case  of  ptl  =  Pi  =  0.5.  In  Fig.  4,  similar  situations  as  in  Fig.  3 
are  considered,  yet  with  p(,  =  0.75  and  p,  =  0.25  at  f  =  0.  The 
solid  line  relates  to  y  =  0  and  follows  Eq.  (20).  Similar  to  Fig.  3, 
the  nonlinear  phase  mismatch  between  the  LG0o  and  LG10 
modes  increases  the  frequency  of  oscillation  of  t  with  respect 
of  L  as  the  relative  nonlinear  coefficient  is  increased.  In  the 
case  of  y  =  3  (dotted),  the  relative  power  transmission  r 


solid,  dashed,  dotted:  y= 0,  0.7,  3 


Fig.  4.  Same  as  Fig.  3,  except  forp0  =  0-75  andpi  =  0.25  at  f  =  0. 


solid,  dashed,  dotted:  7=0,  0.7,  3 


Fig.  5.  Same  as  Fig.  3,  except  forp0  =  0.25  and^q  =  0.75  at  f  =  0. 

remains  nearly  equal  to  unity  for  all  values  of  L.  For  large  val¬ 
ues  of  7  (e.g.,  7  =  3),  the  power  coupling  between  LG0o  and 
LG10  modes  is  inefficient  and  the  amplitude  of  the  oscillation 
in  the  relative  power  carried  by  each  mode  as  a  function  of  L  is 
small  and  is  a  decreasing  function  of  7.  In  other  words,  for  a 
sufficiently  large  values  of  7,  p0  and  p,  remain  unchanged  as 
they  propagate  through  the  GIMF.  The  cumulative  phase  of 
the  LG0o  and  LG10  modes  is  also  nearly  identical  for  large  7 
as  the  two  modes  propagate  along  the  GIMF.  The  combined 
effects  of  the  nearly  unchanged  relative  power  and  the  nearly 
identical  cumulative  phase  of  the  modes  results  in  the  near 
unity  value  of  r  in  this  situation. 

In  Fig.  5,  similar  situations  as  in  Fig.  4  are  considered,  yet 
with  p(,  =  0.25  and  p,  =  0.75  at  C  =  0.  The  solid  line  relates  to 
7  =  0  and  is  identical  to  the  case  of  p„  =  0.75  and  p,  =  0.25 
presented  in  Fig.  4,  in  agreement  with  Eq.  (20).  Similar  to 
Fig.  4,  increasing  the  relative  nonlinear  coefficient  results 
in  a  larger  oscillation  frequency  of  r  with  respect  to  L.  How¬ 
ever,  in  complete  contrast  to  the  behavior  observed  in  Fig.  4 
for  7^0,  increasing  7  initially  lowers  the  minimum  value  of  r 
and  for  sufficiently  large  values  of  7,  the  minimum  can  even  be 
equal  to  zero.  If  7  is  further  increased,  the  minimum  value  of  r 
increases  again  and  eventually  settles  asymptotically  at  a 
value  that  can  be  substantially  different  from  zero.  This  behav¬ 
ior  is  dictated  by  a  very  efficient  exchange  of  power  between 
LG00  and  LG10  modes  as  well  as  considerable  difference  be¬ 
tween  the  cumulative  phases  of  the  modes;  however,  at  a  large 
7  the  cumulative  phases  of  the  modes  become  nearly  identical. 

As  discussed  above,  the  dynamics  of  power  coupling  be¬ 
tween  the  LG0o  and  LG10  modes  is  quite  complex.  In  Fig.  6, 
p„(C  =  L)  is  plotted  as  a  function  of  p0(C  =  0)  for  L  =  Ifwr. 


Fig.  6.  This  figure  shows  the  exchange  of  power  between  the  LG0o 
and  LG10  modes,  where  p0  at  the  output  of  the  GIMF  is  plotted  as  a 
function  of  p0  at  the  input  for  L  =  4.bn. 
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Fig.  7.  Same  as  Fig.  6,  except  for  L  =  50?r,  and  different  choices  of  /. 

In  the  linear  case,  LG0o  and  LG10  modes  are  totally  uncoupled 
and  p0(C  =  L)  =  p0(£  =  0),  as  shown  with  the  diagonal 
dashed  line  in  Fig.  6.  However,  for  y  =  2,  nonlinearity  couples 
the  power  between  LG00  and  LG10  modes  and  p0(c  =  L)  oscil¬ 
lates  as  a  function  of  p0(£  =  0),  shown  with  the  dotted  line  in 
Fig.  6.  It  can  be  seen  thatp0(C  =  L)  oscillates  more  rapidly  as  y 
increases.  In  order  to  isolate  the  nonlinear  effect,  the  disper¬ 
sive  term  rx  is  set  to  zero  and  the  power  coupling  curve  for 
y  =  2  is  re-plotted  as  the  solid  line  in  Fig.  6.  Comparing  the 
solid  and  dotted  lines,  it  can  be  concluded  that  the  effect 
of  the  dispersive  term  )  is  to  tame  the  oscillations, 

especially  for  p0  near  zero  or  unity. 

For  a  longer  GIMF  (at  fixed  y),  the  taming  effect  of  the  dis¬ 
persive  term  will  be  even  more  pronounced  and  the  amplitude 
of  oscillations  will  be  smaller.  In  Fig.  7,  where  L  is  assumed 
to  be  50;r,  the  “artificial”  case  of  c,  =  0,  y  =  0.5  (dotted) 
oscillates  rapidly  and  with  large  amplitude.  However,  for 
the  real  situation  where  r1  «  2  (and  y  =  0.5),  the  rapid  oscil¬ 
lation  is  replaced  with  a  smooth  change  (dashed  line)  that 
does  not  move  far  away  from  the  diagonal  virtual  line  of 
p0(C  =  L)  =  Po(C  =  0).  Once  the  nonlinearity  is  increased  to 
y  =  1.0,  the  oscillations  reappear,  yet  with  a  much  lower  fre¬ 
quency  and  smaller  amplitude  (solid  line)  compared  with  the 
artificial  case  of  /q  =  0.  It  is  noteworthy  that  the  oscillations  in 
the  solid  curve  are  mainly  located  in  the  region  of  the  graph 
where  the  dotted  curve  does  not  oscillate  much. 

In  order  to  examine  the  effect  of  the  total  power  injected  in 
the  GIMF  on  the  relative  power  transmission,  the  transmission 
is  plotted  as  a  function  of  y  for  the  case  of  L  =  4.5 n  in  Fig.  8. 
The  length  of  the  GIMF  is  chosen  in  such  a  way  that  for  the 
linear  case  of  y  =  0,  relative  power  transmission  r  is  at  its 


solid,  dashed,  dotted:  p0= 0.5,  0.25,  0.75 


Fig.  8.  Relative  power  transmission  is  plotted  as  a  function  of  y  (thus 
the  total  power)  for  L  =  4.5;r  for  the  case  of  LGoo  and  LG10  modes, 
whenp0  =  0.5,  0.25,  0.75  at  f  =  0,  in  solid,  dashed,  and  dotted  lines, 
respectively. 


minimum  value  of  rmjll  =  1  -  4p0Pi  for  the  three  cases  of 
p0  =  0.5,  0.25,  0.75  shown  in  Fig.  8.  In  agreement  with  our  pre¬ 
vious  discussions,  increasing  the  power  (increasing  y)  results 
in  an  increase  in  the  value  of  r  and  the  maximum  transmission 
of  rmax  =  1  is  obtained  for  y  k,  0.65,  nearly  independent  of  the 
initial  value  of  p0.  As  y  is  further  increased,  r  decreases  again 
and  follows  an  oscillatory  form  as  a  function  of  y. 

In  practical  devices,  the  GIMF  is  often  considerably  longer 
than  (/?(!)  -  [)( 2»)-1;  therefore,  L  »  1.  In  Fig.  9,  a  scenario 
identical  to  that  explored  in  Fig.  8  is  considered,  except  for 
a  much  longer  GIMF  with  A  =  100.5^.  The  main  difference  be¬ 
tween  the  cases  of  L  =  100.5^- in  Fig.  9andL  =  4.5tr  in  Fig.  8  is 
that  much  lower  power  (much  smaller  y)  is  required  for  the 
GIMF  with  L  =  100.  Fur  to  switch  the  relative  power  transmis¬ 
sion  from  rmm  to  rm;lx .  For  the  case  of  p0  =  0.25,  the  required  y 
for  power  transmission  switching  is  as  low  as  «0.026  as  can  be 
seen  in  the  dashed  line  in  Fig.  9.  Another  important  difference 
between  the  cases  of  L  =  4.5x  and  L  =  100. Sot  is  that  in  the 
latter  case,  the  value  of  y  at  which  rmax  is  obtained  is  very  sen¬ 
sitive  to  the  initial  distribution  of  the  power  among  the  modes. 

B.  Nonlinear  MMI  With  Through  Modes 

In  this  subsection,  the  analysis  of  Subsection  3A  is  extended 
to  the  more  realistic  case  of  five  propagating  modes.  As  men¬ 
tioned  before,  this  is  a  realistic  scenario  for  the  case  of  a 
C-GIMF  at  1550  nm  wavelength.  The  general  observations 
in  this  subsection  and  the  similarities  and  differences  with 
the  two-mode  scenario  of  Subsection  3A  should  cover  the 
overall  nonlinear  dynamics  in  most  realistic  cases  of  NL-MMI 
in  GIMFs. 

In  addition  to  the  difference  between  the  number  of  modes 
in  this  subsection  versus  Subsection  3 .A,  here  the  initial  exci¬ 
tation  amplitudes  are  chosen  such  that  they  satisfy  Eq.  (17) 
(otherwise  the  problem  would  become  intractable).  The  main 
advantage  of  this  choice  is  that  it  complies  with  the  configu¬ 
ration  of  Fig.  1  and  can  be  tested  experimentally.  Another  im¬ 
portant  difference  is  that  in  Subsection  3. A  it  was  assumed 
thatpo  +  P\  =  1,  which  resulted  in  rmax  =  1.  In  a  realistic  case 
where  Eq.  (17)  is  applied  and  the  number  of  guided  modes  is 
finite,  some  of  the  power  coming  from  the  input  SMF  is 
coupled  to  the  radiation  and  cladding  modes.  Therefore,  it 
can  be  shown  [9]  that  the  maximum  transmission  is  given  by 

^max=(l-'P2P)2<l.  (21) 

where  VF  is  defined  in  Eq.  (17)  and  P  is  the  total  number  of 
propagating  modes  with  m  =  0:  P  =  5  in  this  subsection. 


solid,  dashed,  dotted:  p0= 0.5,  0.25,  0.75 


Fig.  9.  Same  as  Fig.  8,  except  for  L  =  100.5^. 
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solid,  dashed,  dotted:  f=0,  0.7,  3 


Fig.  10.  Relative  power  transmission  is  plotted  as  a  function  of  the 
normalized  GIMF  length  L  for  the  case  of  LG00  through  LG40  modes 
(five  zero  angular  modes)  when  p0  =  0-5  at  f  =  0  (>7  =  3  4-  n/8),  for 
y  =  0  (solid),  y  =  0.7  (dashed),  and  y  =  3  (dotted).  The  results  should 
be  compared  with  Fig.  3  where  only  two  modes,  LG0o  and  LG10,  were 
considered. 

Therefore,  in  all  the  plots  shown  in  this  subsection,  rmax  <  1 
and  is  given  by  Eq.  (21). 

In  Fig.  10,  the  relative  power  transmission  is  plotted  as  a 
function  of  the  normalized  GIMF  length  L,  for  //  =  3  +  \/8, 
which  results  in  p0  =  0.5.  Similar  to  the  two-mode  scenario 
plotted  in  Fig.  3,  the  cases  of y  =  0  (solid),  y  =  0.7  (dashed), 
and  y  =  3  (dotted)  are  analyzed.  The  linear  case  (y  =  0)  in 
Fig.  If)  is  very  similar  to  the  two-mode  case  of  Fig.  3,  except 
that  rmln  and  rm;lx  are  different  due  to  the  power  leakage  to  the 
radiation  and  cladding  modes,  as  explained  above.  However, 
the  major  difference  is  that  at  L  =  l.ozr,  y  =  0.7  is  not  nearly 
enough  to  switch  the  relative  power  transmission  from  rrnin  (at 
y  =  0)  to  rmax  in  the  five-mode  scenario.  In  other  words, 
higher  power  is  required  in  the  realistic  five-mode  scenario 
for  transmission  switching  compared  with  the  two-mode  sce¬ 
nario.  This  is  intuitively  expected,  as  the  presence  of  the 
higher  modes  results  in  the  distribution  of  the  intensity  over 
a  larger  cross  sectional  area  of  the  GIMF. 

Next,  the  case  of  p0  =  0.75  for  the  five-mode  scenario  is 
considered  in  Fig.  IT,  which  should  be  compared  with  the 
two-mode  scenario  in  Fig.  4.  According  to  Eq.  (17),  this  case 
can  be  obtained  by  choosing  //  =  3.  Unlike  the  case  of  p0  =  0.5 
studied  above,  there  is  a  strong  similarity  between  the  five¬ 
mode  and  two-mode  scenarios  in  the  p0  =  0.75  case.  At 
L  =  4.5  it,  y  =  0.7  is  nearly  sufficient  to  switch  the  relative 
power  transmission  from  rmm  (at  y  =  0)  to  «rmax.  This  result 


solid,  dashed,  dotted:  y= 0,  0.7,  3 


Fig.  11.  Same  as  Fig.  10,  except  forp0  =  0-75  at  £  =  0  ()/  =  3).  The 
results  should  be  compared  with  Fig.  4  where  only  two  modes,  LG00 
and  LG10,  were  considered. 
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solid,  dashed,  dotted:  y =0,  0.7,  3 


Fig.  12.  Same  as  Fig.  10,  except  for  p0  =  0.25  at  f  =  0 
(i/  =  7  +  \/48).  The  results  should  be  compared  with  Fig.  5  where  only 
two  modes,  LG0o  and  LG10,  were  considered. 

is  quite  important,  because  17  =  3  is  more  readily  achievable 
when  using  standard  commercial  fibers  than  such  large  values 
as  17  =  3  +  \/8  «  5.8.  For  example,  as  mentioned  before,  cou¬ 
pling  C-GIMF  to  SMF-28  at  1550  ran,  results  in  r\  »  2.2,  or  cou¬ 
pling  1060XP  fiber  from  Thorlabs  catalog  to  C-GIMF  leads  to 
11 «  4.27  at  1060  rail  wavelength.  The  major  difference  between 
the  two-mode  and  five-mode  scenarios  is  that  at  the  higher 
value  of  y  =  3,  r  experiences  a  much  higher  amplitude  oscil¬ 
lation  as  a  function  of  L  in  the  five-mode  case.  The  large  am¬ 
plitude  oscillations  in  the  five-mode  scenario  for  y  =  3  are  due 
to  the  efficient  power  transfer  from  the  lowest  order  mode  to 
the  other  four  modes  Each  mode  accumulates  a  different  lin¬ 
ear  and  nonlinear  phase  as  it  propagates  along  the  GIMF  re¬ 
sulting  in  large  variations  in  the  relative  power  transmission, 
which  is  a  sensitive  function  of  both  the  amplitude  and  the 
phase  of  each  mode  [see  Eq.  (15)]. 

Finally,  the  case  of  p0  =  0.25  for  the  five-mode  scenario  is 
considered  in  Fig.  12,  which  should  be  compared  with  the  two¬ 
mode  scenario  in  Fig.  5.  This  case  corresponds  to 
17  =  7  +  \/48  «  13.9,  which  is  very  hard  to  obtain  using  con¬ 
ventional  fibers.  Substantial  power  leakage  into  the  radiation 
and  cladding  modes  is  observed  in  the  five-mode  scenario  and 
the  results  contrast  sharply  from  the  two-mode  scenario  in 
Fig.  5.  Because  of  the  finite  number  of  core-guided  modes  sup¬ 
ported  by  the  finite  GIMF  core  [see  Eq.  (6)],  which  is  trun¬ 
cated  by  the  cladding  at  radius  R,  some  of  the  power  from 
the  input  SMF  does  not  couple  to  the  core-guided  modes 
and  couples  to  the  cladding  and  radiation  modes,  as  estimated 
in  Eq.  (21)  and  also  discussed  in  [9].  The  power  leakage  can  be 
substantial  for  ij  »  1,  as  observed  in  Fig.  12.  There  is  hardly 
any  difference  between  the  linear  case  (y  =  0)  and  y  =  0.7  in 
Fig.  12;  however,  for  y  =  3,  the  relative  power  transmission 
nearly  saturates  over  the  entire  value  of  L  at  rm.lx. 

In  summary,  the  NL-MMI  related  to  p0  =  0.75  in  the  five¬ 
mode  scenario  closely  resembles  that  of  the  two-mode  sce¬ 
nario  presented  in  Subsection  3 .A.  This  similarity  is  lost  as 
the  value  ofp0  is  lowered  (via  increasing  17);  however,  a  larger 
value  of  p0  is  generally  more  accessible  using  commercially 
available  optical  fibers.  Fortunately,  the  power  switching  dy¬ 
namics  (from  rmjn  to  rmax)  is  more  desirable  forp0  =  0-75  than 
p0  =  0.25  in  the  SMF-GIMF-SMF  geometry  shown  in  Fig.  1. 

4.  SATURABLE  ABSORBER  USING  NL-MMI 

In  Section  3,  the  NL-MMI  behavior  of  the  SMF-GIMF-SMF 
geometry  of  Fig.  1  was  discussed  in  great  detail.  Here,  the 
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results  will  be  used  to  explore  the  application  of  this  geometry 
as  an  SA.  Specifically,  it  is  desirable  for  the  SMF-GIMF-SMF 
configuration  to  attenuate  low  power  signals  (y  ~  0),  but  allow 
the  higher  power  signals  to  go  through.  All  analysis  performed 
in  this  section  is  for  the  realistic  case  of  five-modes  (LG0o 
through  LG40),  subject  to  the  constraint  of  Eq.  (17)  for  the 
initial  excitation  amplitudes. 

Figure  13  depicts  the  behavior  of  the  relative  power  trans¬ 
mission  t  as  a  function  of  y  for  a  fixed  value  of  L  =  4.5n  in  a 
SMF-GIMF-SMF  configuration,  for  p0  =  0.5,  0.25,  0.75  at 
f  =  0,  in  solid,  dashed,  and  dotted  lines.  The  value  of  L  is 
chosen  such  that  in  the  linear  case,  the  relative  power  trans¬ 
mission  is  at  its  minimum  value.  In  Subsection  3.B,  it  was 
pointed  out  that  the  case  of  p0  =  0.75  is  the  most  interesting 
scenario  from  an  experimental  point  of  view.  For  p0  =  0.75, 
the  relative  power  transmission  increases  from  rmjll  to  rmax  » 
1  monotonically  as  y  is  increased  from  zero  to  y  ~  0.75.  Be¬ 
yond  y  ~  0.75,  r  goes  through  a  few  low  amplitude  oscillations 
and  saturates  at  rmax  «  1.  This  is  nearly  an  ideal  scenario  for 
an  SA,  where  the  only  downside  with  this  design  is  that  rmhl  is 
considerably  larger  than  zero. 

For  p0  =  0.5  in  Fig.  13,  a  slightly  higher  value  of  y  is  re¬ 
quired  compared  with  p0  =  0.75  to  achieve  the  maximum 
transmission;  moreover,  more  oscillations  are  observed  be¬ 
yond  the  transmission  peak  value,  which  may  not  be  a  desir¬ 
able  feature  for  an  SA.  Another  important  feature  of  this  plot  is 
the  large  low-transmission  plateau  and  the  sudden  rise  of  z 
near  the  peak  values,  which  is  desirable  for  pulse  shortening 
in  an  SA.  The  case  of  p0  =  0.25  is  also  shown  in  Fig.  13;  the 
peak  value  of  transmission  in  this  scenario  is  only  near  50% 
due  to  power  coupling  to  the  radiation  and  cladding  modes, 
as  discussed  in  Subsection  3T5,  making  this  an  undesirable  de¬ 
sign.  The  results  in  Fig.  hi  should  be  compared  with  the  two- 
rnode  scenario  in  Fig.  8.  A  notable  feature  in  Fig.  8  is  the  large- 
amplitude  oscillations  of  r  as  a  function  of  y;  this  behavior  is 
very  different  from  Fig.  13,  where  z  was  nearly  saturated  to 
Tmax  beyond  y  a  0.75  for  p0  =  0.75. 

In  Fig.  M,  a  similar  case  to  the  one  discussed  in  Fig.  IB  is 
considered,  except  for  a  longer  GIMF  section  where 
L  =  100.5®.  Similarly,  this  value  of  L  is  chosen  such  that,  in 
the  linear  case,  the  relative  power  transmission  is  at  its  mini¬ 
mum  value.  The  most  notable  differences  between  the  case  of 
L  =  100.5®  in  Fig.  14  and  L  =  4.5 ®  in  Fig.  13  are  that  rmax  is 
obtained  at  a  lower  value  of  y  for  the  longer  GIMF  design.  This 
is  somewhat  expected,  given  that  the  nonlinear  phase  is 


solid,  dashed,  dotted:  p0=0. 5,  0.25,  0.75 


Fig.  13.  Relative  power  transmission  is  plotted  as  a  function  of  y 
(thus  the  total  power)  for  L  =  4.5®  for  the  case  of  five  modes, 
whenp0  =  0.5,  0.25,  0.75  at  f  =  0,  in  solid,  dashed,  and  dotted  lines, 
respectively. 


solid,  dashed,  dotted:  p0= 0.5,  0.25,  0.75 


Fig.  14.  Same  as  Fig.  13,  except  for  L  =  100.5®. 

cumulative.  However,  comparing  the  five-mode  scenario  in 
Fig.  14  with  the  two-mode  case  of  Fig.  9  shows  that  a  much 
higher  value  of  y  is  required  in  the  five-mode  scenario  to  ob¬ 
tain  rmax;  therefore,  the  five-mode  and  two-mode  scenarios 
differ  substantially  from  each  other. 

In  Figs.  13  and  14,  the  value  of  L  is  chosen  such  that  z  =  rlllm 
is  obtained  for  y  =  0.  While  this  may  seem  a  reasonable 
choice,  this  is  not  required  from  an  SA;  rather,  all  that  is  re¬ 
quired  from  an  SA  is  to  efficiently  discriminate  against  low 
power  signals  with  minimal  impact  on  the  high  power.  There¬ 
fore,  it  is  necessary  to  examine  other  possibilities  for  L  and 
find  the  most  efficient  SA  that  is  possible  with  the  SMF- 
GIMF-SMF  geometry  of  Fig.  1.  A  relevant  metric  for  the  pulse 
power  discrimination  in  an  SA  is  the  value  of  ( dz/dy ),  which 
characterizes  the  sensitivity  of  an  SA  to  a  change  in  the 
pulse  power. 

Figure  15  shows  (5z/z),  defined  as 

(<5r/r)  =  -(r|j;=o.ooi  ~  Tlp=o)'  (22) 

z 

as  a  fimction  of  the  normalized  GIMF  length.  (f>z/z)  is  used  as 
the  metric  because  it  is  related  to  {dz/dy)  via  (5z/z) « 
z~l(dz  /  dy)8y  for  8y  =  0.001.  The  normalized  length  of  the 
GIMF  is  chosen  to  be  L  =  300®  +  8L,  and  the  horizontal  axis 
in  Fig.  Hi  is  expressed  as  8L  in  the  range  [— zr/2,  +®/2].  The 
reason  for  this  limited  range  in  8L  is  that  (8r/z)  repeats  in 
a  nearly  periodic  fashion  outside  this  range,  unless  SL  be¬ 
comes  very  large.  It  will  be  discussed  later  in  Section  5  that 


Fig.  15.  Relative  change  in  the  normalized  power  transmission  in  the 
SMF-GIMF-SMF  geometry  defined  in  Eq.  (22)  plotted  as  a  function  of 
SL,  where  the  length  of  the  GIMF  section  is  given  by  L  =  300®  +  SL.  A 
large  positive  value  of  ( Sz/z )  is  desirable  for  saturable  absorption.  The 
plots  are  presented  for  r\  =  2  (solid)  and  r\  =  3  (dashed). 
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Fig.  16.  Relative  power  transmission  is  shown  as_  function  of  SL , 
where  the  length  of  the  GIMF  section  is  given  by  L  =  300#  +  SL.  r 
presented  in  this  figure  is  the  average  value  between  r|^=o.ooi  and 
t\y=o-  This  figure  is  intended  to  show  that  the  value  of  r  is  not  too 
low  near  the  peak  value  of  ( Sz/z )  in  Fig.  15.  The  plots  are  presented 
for  //  =  2  (solid)  and  )/  =  3  (dashed). 

even  a  seemingly  small  value  of  y  =  0.001  can  translate  to 
multi-kilowatts  of  peak  power  in  some  commercial  GIMFs, 
such  as  C-GIMF.  Therefore,  our  main  intention  in  choosing 
Sy  =  0.001  in  Eq.  (22)  and  subsequent  figures  is  to  ensure  that 
our  design  and  observations  are  relevant  for  practical 
situations. 

For  both  17  =  2  (solid  line)  and  r\  =  3  (dashed  line)  in  Fig.  J_5, 
the  maximum  value  of  (Sz/z)  is  obtained  at  SL  «  -0.4,  while 
the  large  low-transmission  plateau  already  observed  in  Figs.  13 
and  14  at  near  y  =  0  is  responsible  for  the  small  value  of  (St /t) 
near  SL  =  0.  The  maximum  value  of  (St /t)  is  around  1%  for 
rj  =  2  and  5%  for  ;/  =  3.  It  should  be  pointed  out  that  SAs 
can  be  used  to  mode-lock  lasers  even  when  their  modulation 
depth  (MD)  is  as  low  as  0.5%  [27].  MD  is  defined  as  the  maxi¬ 
mum  change  in  transmission  in  an  SA.  Although  MD  is  often 
the  standard  quantity  used  to  assess  the  performance  of  SAs, 
the  large  values  of  MD  obtained  in  the  SMF-GIMF-SMF 
(between  rmill  and  rm;lx )  can  be  misleading,  because  a  very 
high  peak  power  may  be  needed  to  access  the  full  range  of 
the  MD  (some  numerical  values  will  be  presented  in  Section  5). 
Therefore,  (St /t)  has  been  used  in  this  paper  and  is  identical 
to  MD  within  the  reasonably  accessible  range  of  pulse  peak 
powers. 

It  is  required  to  ensure  that  the  total  transmission  through 
the  SMF-GIMF-SMF  geometry  is  not  too  low;  in  other  words, 
the  SA  is  not  too  lossy  near  the  peak  value  of  (5t/t).  In  Fig.  16, 
the  average  value  (t^o.ooi  +  'Z‘lp=o) /2  is  plotted  as  a  function 
of  SL  for  77  =  2  (solid  line)  and  //  =  3  (dashed  line).  It  can  be 


L=300#+SL 


seen  that  the  average  relative  power  transmission  is  reason¬ 
ably  large  near  the  peak  value  of  (St/ t)  of  Fig.  15;  therefore, 
the  SA  can  be  regarded  as  a  viable  design. 

In  Fig.  F7,  (Sz/z)  is  plotted  in  the  same  fashion  as  in  Fig.  15, 
except  for  r\  =  4.  The  reason  separate  figures  are  used  to  plot 
11  =  2,  3  and  17  =  4  is  that  the  vertical  scales  are  different.  As 
can  be  seen  in  Fig.  17,  (St /t)  as  large  as  15%  is  possible  in  this 
case.  While  this  may  seem  like  a  highly  desirable  SA,  the  rel¬ 
ative  power  transmission  plot  for  1/  =  4  in  Fig.  18  shows  that  r 
is  quite  small  near  the  peak  values  of  (Sz/z).  Perhaps  an 
optimum  operation  of  this  device  can  be  accomplished  at 
SL  «  -0.3  where  5t/t  «  9%. 

5.  DISCUSSION  AND  CONCLUSION 

In  this  paper,  a  detailed  analysis  of  the  NL-MMI  behavior  of  the 
SMF-GIMF-SMF  geometry  shown  in  Fig.  1_  has  been  presented. 
Recently,  this  setup  was  successfully  used  to  create  very  low- 
loss  couplers  between  two  SMFs  with  very  different  mode-field 
diameters  [9,10].  Here,  it  has  been  shown  that  the  transmission 
through  this  coupler  has  a  very  interesting  nonlinear  behavior 
and  can  be  potentially  used  for  switching  purposes  or  for  satu¬ 
rable  absorption  in  all-fiber  mode-locked  lasers. 

The  key  parameters  that  affect  the  NL-MMI  in  this  geometry 
are  the  ratio  of  the  mode-field  diameter  of  the  LG0o  mode  in 
GIMF  to  the  mode-field  diameters  of  the  SMFs  [characterized 
by  17  defined  in  Eq.  (17)];  the  length  of  the  GIMF  section  (char¬ 
acterized  by  the  relative  length  parameter  L);  the  relative  non¬ 
linear  coefficient  y,  which  incorporates  the  total  power  and 
nonlinear  coefficient  of  the  GIMF;  and  the  total  number  of 
propagating  modes  in  the  GIMF. 

It  was  observed  that  the  best  performances  for  switching 
and  intensity  discrimination  applications,  such  as  in  SAs, 
can  be  obtained  when  the  value  of  17  defined  in  Eq.  (17)  is 
not  too  large;  17  <  4  seems  like  a  reasonable  choice.  This  is 
easily  achievable  for  reasonable  selection  of  commercially 
available  fibers. 

Like  most  other  nonlinear  devices,  the  interaction  length 
plays  a  very  important  role  in  the  behavior  of  the  SMF- 
GIMF-SMF  device.  The  dimensionless  normalized  length  L 
has  been  used  throughout  this  paper,  which  is  the  length  of 
the  GIMF  nonnalized  by  /?(1)  -  /?(2)-  This  normalization  factor 
is  (183  71m)  1  for  a  C-GIMF  operating  at  1550  ran.  Therefore, 
L  =  4.5#,  L  =  100.5#,  and  L  =  300#  used  in  this  paper  trans¬ 
late  to  L  =  2.6  mm,  5.78,  and  17.24  cm  for  the  length  of  the 
GIMF  section,  which  are  all  reasonable  lengths  for  device  ap¬ 
plications.  For  SA  purposes,  L  =  300#  +  SL  was  explored  in 
Figs.  15  and  17,  where  it  was  realized  that  some  level  of  fine 
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tuning  of  the  order  of  SL  ~  n  is  required,  which  is  SL  ~  575  pm 
for  C-GIMF  at  1550  nm  wavelength.  This  situation  is  not  differ¬ 
ent  from  that  of  Refs.  [9, 10],  and  the  fact  that  the  required  fine- 
tuning  range  is  so  small  means  that  it  can  be  easily  achieved 
by  polishing  the  GIMF  to  the  desired  length.  It  should  be  em¬ 
phasized  that  the  problem  of  fine-tuning  is  not  very  serious 
and  the  length  needs  to  be  adjusted  over  a  maximum  one 
period  of  SL  ~  n,  because  reducing  LtoL-n  shifts  the  peri¬ 
odic  transmission  pattern  by  one  full  period.  Other  pragmatic 
approaches  can  be  taken  instead  of  fme-timing  the  length.  For 
example,  an  ensemble  of  twenty  SMF-GIMF-SMF  couplers  can 
be  built  around  the  same  average  GIMF  length.  Using  this  ap¬ 
proach,  the  value  of  SL  will  be  uniformly  distributed  over  the 
range  of  [-n/2,  +n/2]  and  one  of  the  elements  of  the  ensemble 
will  likely  perform  up  to  the  required  specification. 

The  relative  nonlinear  coefficient  y  is  defined  in  Eq.  (12). 
For  C-GIMF  at  1550  nm  wavelength,  y  =  10  corresponds  to 
a  power  of  P  «  10.0  MW,  given  that  the  nonlinear  coefficient 
y  (defined  for  the  LG0o  mode  in  GIMF)  is  smaller  than  that  of 
the  SMF-28  by  a  factor  // «  2.2.  For  practical  fiber-based  mode¬ 
locking  applications,  this  value  is  unacceptably  large.  The  fea¬ 
sible  range  of  operation  for  the  peak  power  of  a  mode-locked 
fiber  laser  is  from  several  to  a  likely  maximum  of  a  few  hun¬ 
dred  KW  of  peak  power.  If  the  GIMF  is  made  from  highly  non¬ 
linear  material,  e.g.,  chalcogenide  glass  with  a  nonlinear 
coefficient  of  up  to  1000  times  larger  than  silica,  then  y  = 
1.0  will  correspond  to  P  «  10.0  KW 

The  feasibility  of  using  the  SMF-GIMF-SMF  geometry  has 
been  shown  for  SA  applications,  even  for  y  =  0.001,  in  Fig.  17, 
which  corresponds  toP«  10.0  KW  using  C-GIMF  and  SMF-28 
commercial  fibers  at  1550  nm  wavelength.  Therefore,  it  is 
probable  that  the  SMF-GIMF-SMF  geometry  can  be  used  even 
with  conventional  commercially  available  fibers  to  mode-lock 
fiber  lasers  at  the  presently  achievable  power  levels.  The  main 
advantage  of  this  geometry  is  that  it  can  be  designed  to  oper¬ 
ate  at  a  much  higher  power  level  compared  with  the  existing 
SAs  that  exploit  nonlinear  polarization  rotation,  semiconduc¬ 
tors,  or  carbon  nanotubes.  Therefore,  it  might  offer  a  very  at¬ 
tractive  solution  for  scaling  up  the  pulse  energy  and  peak 
power  in  mode-locked  fiber  lasers  that  are  currently  limited 
by  the  existing  SA  technology. 

Throughout  this  paper,  it  was  observed  that  having  a 
smaller  number  of  propagating  modes  in  the  GIMF  is  benefi¬ 
cial,  especially  in  the  design  of  an  SA.  While  specific  designs 
need  to  be  carefully  analyzed  using  the  formalism  laid  out  in 
this  manuscript,  this  can  be  taken  as  a  general  guideline  and 
highly  multimode  GIMFs  should  be  avoided  if  possible.  We  do 
not  think  that  the  quality  and  the  telecommunication  band¬ 
width  of  the  GIMF  (e.g.,  due  to  the  centerline  defect  [28]) 
has  any  appreciable  effect  on  the  NL-MMI  behavior  and  the 
SA  performance  of  the  SMF-GIMF-SMF  geometry. 

Finally,  we  would  like  to  comment  on  the  possibility  of  ex¬ 
citing  nonzero  angular  momentum  modes,  which  are  excluded 
from  our  analysis.  It  is  possible  to  excite  m  =  ±  1  modes  via 
coupling  them  to  m  =  0  modes  by  bending  the  GIMF.  In  order 
for  the  power  exchange  to  operate  efficiently  in  a  resonant 
fashion,  the  phase  mismatch  between  the  m  =  0  and  m  = 
±1  must  be  compensated  for  by  the  bend;  e.g.,  long-period 
gratings  are  optimized  to  provide  this  phase-matching  condi¬ 
tion  [17],  However,  in  a  GIMF,  the  m  =  0  and  m  =  ±1  modes 
always  belong  to  different  mode  groups  with  a  very  large 


mismatch  in  their  propagation  constants,  according  to 
Eq.  (4).  Therefore,  a  very  small  bending  radius  (submillimeter 
for  C-GIMF)  is  required  for  efficient  phase  matching.  There¬ 
fore,  coupling  between  m  =  0  and  m  =  ±  1  will  not  be  an  ef¬ 
ficient  process  in  GIMFs,  unless  assisted  by  specialty  devices 
such  as  long-period  gratings  [17], 

Coupling  between  between  m  =  0  and  m  =  ±2  is  induced 
via  the  random  variations  in  the  diameter  of  the  GIMF,  caused 
by  manufacturing  errors  [29],  and  can  even  happen  within  the 
same  mode  group  with  minimal  phase  matching  issue.  While 
this  can  be  minimized  by  selecting  a  high  quality  GIMF,  in 
practice  this  should  not  cause  any  issues  for  GIMFs  shorter 
than  approximately  one  meter  in  length. 

It  has  also  been  previously  shown  that  fabrication  errors 
due  to  the  misalignment  of  the  fibers  can  be  minimized  to 
the  extent  where  they  do  not  have  a  noticeable  impact  on 
the  linear  operation  of  the  device  [10],  Similar  robustness  is 
expected  in  the  nonlinear  regime  of  operation. 
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Abstract:  The  main  differences  in  nonlinear  switching  behavior  between 
multicore  versus  multimode  waveguide  couplers  are  highlighted.  By 
gradually  decreasing  the  separation  between  the  two  cores  of  a  dual-core 
waveguide  and  interpolating  from  a  multicore  to  a  multimode  scenario,  the 
role  of  the  linear  coupling,  self-phase  modulation,  cross-phase  modulation, 
and  four-wave  mixing  terms  are  explored,  and  the  key  reasons  are  identified 
behind  higher  switching  power  requirements  and  lower  switching  quality  in 
multimode  nonlinear  couplers. 
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1.  Introduction 

The  nonlinear  switching  (NLS)  properties  of  waveguide  arrays  and  multicore  fiber  couplers 
have  recently  attracted  considerable  attention  for  mode-locked  fiber  laser  applications  [1-5]. 
In  the  linear  regime,  where  optical  power  is  low,  neighboring  waveguides  exchange  optical 
power  periodically;  the  linear  coupling  is  caused  by  the  modal  overlap  of  adjacent  waveguides 
and  is  most  efficient  when  adjacent  modes  have  identical  propagation  constants  [6].  In  the 
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nonlinear  regime,  where  optical  power  is  high,  nonlinear  effects  alter  the  refractive  index  of 
each  waveguide  and  consequently  detune  the  effective  propagation  constants  of  the  modes, 
reducing  the  power  exchange  efficiency  between  adjacent  cores  [7].  By  only  retaining  the  light 
in  the  launch  waveguide  at  the  output,  it  is  possible  to  achieve  power-dependent  transmission 
and  intensity  discrimination  required  for  mode-locked  operation  [1,2]. 

Recently,  it  was  shown  that  nonlinear  multimodal  interference  in  a  graded-index  multimode 
fiber  has  intensity  discrimination  properties  as  well  and  can  be  used  for  NLS  [8].  However,  it 
was  observed  that  a  much  higher  power  is  required  for  NLS  based  on  the  nonlinear  multi-modal 
interference  in  multimode  fibers  compared  to  the  NLS  based  on  the  nonlinear  mode  coupling 
in  multicore  fibers.  Moreover,  the  NLS  quality  of  the  multimode  fiber  was  shown  to  be  inferior 
to  those  reported  for  multicore  fibers. 

The  objective  of  this  manuscript  is  to  highlight  the  differences  between  the  NLS  behavior  in 
multimode  versus  multicore  waveguides.  Several  key  questions  are  addressed:  Why  is  higher 
power  required  to  achieve  NLS  in  multimode  versus  multicore  waveguides?  What  is  the  main 
reason  behind  the  inferior  NLS  quality  in  multimode  waveguides?  Is  it  because  of  the  cross¬ 
phase  modulation  (XPM)  and  four- wave  mixing  (FWM)  terms?  Or  because  of  the  way  the 
modes  are  excited  in  the  multimodal  junctions? 


Fig.  1.  (a)  The  refractive  index  profile  of  the  double-waveguide  nonlinear  coupler  is  shown. 

(b)  For  weakly  coupled  waveguides  ( d  =  10A),  the  even  and  odd  mode  profiles  Fe(x)  and 
F0(x)  in  Pictures#  are  sketched  above  the  index  profile  and  F\  ( x )  and  Fijx)  in  Picture- SB 
are  sketched  below  the  index  profile,  (c)  is  similar  to  (b),  except  the  waveguides  are  strongly 
coupled  for  d  =  2\.  Results  are  shown  for  the  TE  polarization,  where  the  transverse  electric 
field  vector  is  pointing  in  the  vertical  direction  in  these  figures. 

In  order  to  present  the  arguments  in  a  concrete  model,  NLS  is  explored  for  the  Transverse 
Electric  (TE)  polarization  in  a  pair  of  identical  one-dimensional  waveguides  shown  in  Fig.  1(a). 
The  design  parameters  are  such  that  each  waveguide  supports  only  a  single  TE-polarization 
mode,  with  nci  =  1.5,  nco  =  1.506,  and  a  =  3/1,  where  A  is  the  optical  wavelength.  We  refer  to 
the  spatial  transverse  mode  profile  of  a  “single”  waveguide  as  Fw(x)  throughout  this  paper. 

In  the  linear  regime,  light  propagation  through  this  double-waveguide  is  commonly  treated 
using  the  standard  coupled-mode  theory,  which  is  valid  as  long  as  the  waveguides  are  weakly 
coupled  [6].  In  Fig.  1(b)  where  the  separation  between  the  waveguides  is  d  =  10A,  the  stan¬ 
dard  coupled  mode  theory  can  be  reliably  applied  to  the  individual  modes  of  the  waveguides 
whose  profiles  are  identified  as  F\  {x)  and  Fjix).  The  overlap  between  these  modal  profiles  de¬ 
termines  their  coupling.  However,  a  rigorous  approach,  which  is  applicable  to  both  weak  and 
strong  coupling,  is  based  on  directly  solving  for  the  even  and  odd  supermodes  of  the  double¬ 
waveguide  and  treating  the  light  propagation  as  a  modal  interference  problem.  The  even  and 
odd  eigenmodes  of  the  full  index  profile  are  sketched  as  Fe(x)  and  F„(x)  in  Fig.  1(b).  When 
the  waveguides  are  only  weakly  coupled,  the  two  approaches  can  be  identically  mapped  to  one 
another  by  F\2  =  (Fe  ±  F„) / \/2.  We  will  take  this  relationship  as  the  definition  of  F\  (x)  and 
F)  (x) .  Therefore,  for  the  case  of  the  weakly  coupled  waveguides,  we  can  write  F\  =  F2  =  Fw, 
where  =  implies  similarity  up  to  a  shift  in  the  x-coordinate. 

The  case  of  the  strongly  coupled  waveguides  with  d  =  2A  is  shown  Fig.  1(c),  where  Fe(x) 
and  Fa(x)  profiles  are  sketched  as  the  exact  eigenmodes  of  the  system.  The  standard  coupled 
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mode  theory  is  not  applicable  here  ,  and  F\  (x)  and  f)  (x)  profiles  obtained  from  F\2  (x)  =  (Fe  ± 
Fn) /\/l  are  far  from  the  eigenmodes  of  the  individual  waveguides.  In  this  case,  F\  ^  f  2  ^  Fw. 

In  this  paper,  the  NLS  behavior  of  the  double-waveguide  of  Fig.  1(a)  is  studied  as  a  function 
of  the  waveguide  separation.  This  analysis  provides  an  interpolation  from  the  multicore  to  the 
multimode  setup  as  the  two  waveguides  are  brought  closer  together  and  eventually  merged. 
The  NLS  problem  is  initially  treated  using  the  generalized  nonlinear  Schrodinger  equation 
(GNLSE)  applied  to  the  even  and  odd  supermodes  of  the  double-waveguide  (Fe  and  F„),  which 
are  the  exact  eigenmodes  of  the  Helmholtz  equation.  This  treatment  is  rigorous  and  is  referred 
to  as  Picture -sX .  The  nonlinear  propagation  problem  is  then  transformed  to  the  language  of  the 
F\  (x)  and  /^(x)  profiles,  referred  to  as  Picture-38.  Although,  /q(x)  and  /^(x)  are  not  eigen¬ 
modes  of  the  individual  waveguides  and  take  their  meaning  only  from  Fe(x)  and  F0 (x),  the 
nonlinear  propagation  problem  can  be  recast  in  their  language  with  no  loss  of  generality  and 
both  Pictures  are  equally  valid.  However,  it  will  be  seen  that  Picture-38  is  more  suitable  for 
the  purposes  of  this  paper,  as  the  role  of  the  XPM/FWM  versus  self-phase  modulation  (SPM) 
terms  and  the  initial  modal  excitations  are  more  clear  in  this  Picture.  This  holds  especially  true 
in  the  limit  of  weakly  coupled  waveguides,  where  the  XPM/FWM  terms  are  completely  absent 
in  the  nonlinear  propagation  equation  in  Picture-38 . 


2.  The  formulation  and  results 


In  what  follows,  only  the  continuous  wave  (CW)  limit  is  considered,  in  order  to  reduce  the 
complexity  of  the  analysis  and  to  ensure  that  the  physics  is  not  buried  under  phenomena  that 
are  not  essential  in  conveying  the  intended  message.  However,  the  temporal  effects  can  be 
easily  included  and  the  main  observation  will  not  be  affected  .  The  observations  are  general  and 
equally  apply,  at  least  qualitatively,  to  Transverse  Magnetic  (TM)  waves,  and  optical  fibers. 

In  the  coupled  waveguide  problem  of  Fig.  1(a),  the  lowest  order  modes  are  a  pair  of  z- 
propagating  bound  states  with  the  even  and  odd  spatial  profiles  Fe(x)  =  Fe(—x)  and  F0(x)  = 
— F0(—x )  being  invariant  along  the  y-axis.  The  symmetries  of  these  two  profiles  are  dictated  by 
the  parity  symmetry  of  the  refractive  index  profile  «(x)  =  n(—x)  and  its  invariance  along  the 
y-axis.  The  GNLSE  describing  the  CW  evolution  of  light  in  the  two-mode  double-waveguide 
in  Picture- sX  can  be  written  as  [9, 10] 


+  i 


'(*?) 


s 

V,fC,p=< 


fnvKpAvAKAp 


H=e,o, 


(1) 


where  /\u  (z)  is  the  slowly  varying  envelope  of  the  electric  field  of  the  pth  mode  (even  or  odd 
mode)  with  the  propagation  constant  fip.  Pp  =  \Ap \ 2  is  the  linear  power  density  carried  by  this 
mode  per  unit  v- length.  The  indices  can  take  the  value  of  e  or  o  corresponding  to  even  and 
odd  modes.  We  also  define  5j8f  =  —8[ 30  =  jie  —  /3ref,  where  j3ref  =  (/3C  +  /30)/2.  The  nonlinear 
coupling  coefficients  are  given  by 

(2) 


fpvKp  —  J  FpFvFKFp  dx , 


where  the  mutually  orthogonal  spatial  profiles  are  assumed  to  be  normalized  according  to 
f  Fjdx  =  f  Fgdx  =  1 .  «2  is  the  nonlinear  index  coefficient,  and  (Oq  is  the  carrier  frequency. 
fpvKp  is  a  fully  symmetric  tensor,  and  feooo  and  feeeo  vanish  due  to  the  parity  symmetry. 

In  Fig.  2(a),  8  jir  is  plotted  as  a  function  of  the  normalized  separation  d/X  between  the 
waveguides.  At  a  large  separation,  the  even  and  odd  modes  are  nearly  degenerate  in  the  value 
of  the  propagation  constant.  However,  8 fir  rapidly  increases  as  the  waveguides  are  brought 
closer  together  and  the  degeneracy  between  the  even  and  odd  modes  is  strongly  broken. 

In  Fig.  2(b),  the  non-zero  nonlinear  coupling  coefficients  of  the  even  and  odd  modes  in 
Picture -sX  are  plotted  as  a  function  of  the  normalized  waveguide  separation.  The  nonlinear 
coupling  coefficients  feeee,feeoo  and  f000o  are  equal  when  the  two  waveguides  are  far  apart  and 
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Fig.  2.  (a)  The  splitting  between  the  propagation  constants  of  the  modes  increases  signif¬ 
icantly  as  the  waveguides  are  brought  closer  together,  (b)  The  nonlinear  coupling  coeffi¬ 
cients  for  the  even  and  odd  modes  are  shown  as  a  function  of  the  normalized  separation, 
where  their  degeneracy  is  removed  when  the  waveguides  are  strongly  interacting. 

hence,  weakly  coupled.  However,  as  the  two  waveguides  get  closer  ,  the  difference  between  the 
nonlinear  couplings  increases.  Other  elements  of  the  fpvKp  tensor  not  shown  in  Fig.  2(b)  are 
related  to  the  plotted  elements  by  symmetry  relationships  of  the  tensor;  e.g.,  feoeo  =  feeoo- 
The  large  values  of  the  XPM/FWM  terms  in  Picture -stf  complicate  the  analysis  of  the  NLS. 
In  Picture -PS,  the  nonlinear  propagation  is  more  intuitive  and  the  XPM/FWM  terms  are  absent 
in  the  limit  of  weakly  coupled  waveguides.  In  Picture -PS,  the  field  amplitudes  are  defined  as 
A]  2  :=  (Ae  ±A0)/y/2.  The  GNLSE  Eq.  (1)  is  transformed  to 

^  =  id  pe A  jt  +  i  Y,  fjklmAkAiA*m,  j  =  1,/  =  2,  or  j  =  2,/  =  1,  (3) 

dz  \  c  /  kjm=12 

where  the  nonlinear  coefficients  fjki,n  are  given  by 

f\  I  1  I  ~  f2222  =  i  feeee  T  foooo  T  ^2 feeoo )  /4 ,  (4) 

fll22  =  ifeeee  A  foooo  ^  feeoo)/  4,  ,/l  I  ]  2  ~  f\  222  —  ifeeee  foooo) /4. 

The  results  of  Eq.  (4)  are  consistent  with  the  field  profiles  defined  in  Picture- SS  as  F\2  = 
(Fe±F0)/V 2,  given  fjki,n  =  f  FjFkFiFm  dx  with  j,k,l,m  =  1,2.  fjUm  is  a  fully  symmetric 
tensor.  It  remains  invariant  if  all  indices  valued  at  1  are  changed  to  2  and  vice  versa;  e.g., 
/l  1 12  =  /1222- 


Fig.  3.  The  nonlinear  coupling  coefficients  in  Picture-SS  are  shown  as  a  function  of  d/X. 
The  SPM  coefficients  in  (a)  are  much  stronger  than  the  XPM/FWM  coefficients  in  (b). 


The  non-zero  elements  of  fjkim  are  plotted  in  Figs.  3(a)  and  3(b).  The  SPM  nonlinear  cou¬ 
pling  terms  fun  and  f2222  are  larger  than  the  XPM/FWM  terms  by  more  than  an  order  of 
magnitude,  even  when  the  two  waveguides  are  merged.  This  behavior  is  dictated  by  the  rela¬ 
tionships  presented  in  Eq.  (4)  and  the  near  degeneracy  of  fpvKp  in  Fig.  2(b),  especially  at  weak 
coupling.  Thus,  in  Picture-PS,  the  nonlinear  switching  is  dominated  by  SPM,  which  is  easier  to 
analyze. 

The  NLS  behavior  of  the  waveguide  coupler  is  shown  in  Fig.  4.  A  multicore  NLS  configura¬ 
tion  is  sketched  in  Fig.  4(a),  where  the  beam  is  injected  into  the  two-core  switch  using  an  input 
coupler  waveguide  and  is  collected  at  an  output  coupler.  The  length  of  the  nonlinear  switch  is 
equal  to  the  half-beat-length  L/„  so  that  at  low  powers,  all  the  injected  light  is  transferred  to  the 
second  waveguide  and  no  power  is  collected  at  the  output.  However,  when  the  input  power  is 
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Fig.  4.  (a)  A  sketch  of  the  NLS  device,  (b)  The  relative  power  transmission  is  plotted 
as  a  function  of  y.  The  solid  (black)  curve  corresponds  to  weakly  coupled  waveguides 
with  d/X  =  10.  The  dashed  (red)  and  dotted  (blue)  NLS  curves  correspond  to  the  strongly 
interacting  merged  waveguide  with  d/X  =  0,  where  in  case-1  the  injected  and  collected 
beam  profiles  are  F\  (x),  while  in  case-2,  the  injected  and  collected  beam  profiles  are  Fw(x). 


increased,  the  nonlinearity  detunes  the  coupling  between  the  two  waveguides  and  all  the  power 
is  collected  at  the  output  port.  This  behavior  is  shown  in  Fig.  4(b),  where  the  transmission  T  is 
plotted  as  a  function  of  the  effective  nonlinear  coefficient  defined  as  y  =  rjCOqP  / ( c8fie ),  where 
P  is  the  linear  optical  power  density  injected  from  the  input  coupler. 

In  Fig.  4(b),  the  solid  (black)  line  shows  the  NLS  behavior  of  the  double-waveguide  coupler 
when  the  normalized  separation  is  d/X  =  10.  In  this  case,  the  injected  and  collected  field  pro¬ 
files  are  assumed  to  be  those  of  a  single  waveguide  identical  to  the  waveguides  in  the  nonlinear 
coupler.  The  relative  transmitted  power  can  be  obtained  from  [8] 


T  = 


1 

P2 


X  a*(0)a;(lc)|2  =  -^|  X  a;(0)a,(lc)\2 


(5) 


=  1,2  *  H=e,o 

where  Lc  is  the  total  length  of  the  nonlinear  switch.  In  this  case,  the  input  boundary  conditions 
(BC)  for  the  GNLSE  of  Eq.  (3)  are  A  j  (0)  / VP  =  1  and  A2  (0)  / VP  =  0.  These  BCs  hold  because 
the  injected  and  collected  beam  profiles  are  Fw(x)  and  is  mode-matched  to  F)  (x)  of  Picture- 33, 
which  is  due  to  the  weak  coupling  between  the  two  waveguides. 

When  the  waveguides  are  brought  closer  together,  F\  (x)  can  no  longer  be  approximated  as 
F„{x),  as  can  also  be  seen  in  Fig.  1(c).  In  Fig.  4(b),  the  dashed  (red)  line  shows  the  NLS  behavior 
of  the  double-waveguide  coupler,  when  the  normalized  separation  is  d/X  =  0,  which  is  the  limit 
of  a  standard  single-core  multimode  coupler.  In  this  case  (case-1),  the  BCs  are  assumed  to  be 
A i  (0) /VP  =  1  and  A2(0) / VP  =  0.  However,  given  the  difference  between  F\  (x)  and  Fw(x)  in 
this  strong  coupling  regime,  it  is  clear  that  these  BCs  cannot  be  achieved  using  an  input  coupler 
waveguide  identical  to  those  used  in  the  nonlinear  switch.  In  fact,  the  injected  and  collected 
beam  profiles  must  be  modified  by  additional  optics  (such  as  a  spatial  light  modulator)  or  must 
come  from  special  waveguides  in  order  to  mode-match  to  F\  (x).  Even  if  Fj(x)  is  used  for  the 
injection,  it  might  be  more  desirable  to  collect  the  beam  directly  from  the  nonlinear  coupler  in 
the  form  of  Fw(x).  In  this  case,  it  is  possible  to  slightly  change  the  device  length  away  from  L/, 
to  improve  the  coupling  efficiency  at  high  power  at  the  expense  of  a  perfect  attenuation  at  low 
power  (linear  case).  This  is  a  reasonable  trade-off  because  efficient  transmission  at  high  power 
is  desired  for  mode-locked  lasers  while  some  level  of  compromise  in  nonlinear  modulation 
depth  is  acceptable. 

The  difference  between  the  dashed  (red)  line  of  d/X  =  0  compared  with  the  solid  (black) 
line  of  d/X  =  10  is  caused  by  the  XPM/FWM  terms;  if  the  XPM/FWM  terms  in  Fig.  3(b)  are 
artificially  set  to  zero,  the  dashed  (red)  line  in  Fig.  4(b)  falls  almost  exactly  on  the  solid  (black) 
line.  The  XPM/FWM  terms  nearly  double  the  switching  threshold  in  units  of  y.  However,  given 
that  Spe  is  nearly  555  times  larger  for  d/X  =0  than  d/X  =  10  [see  Fig.  2(a)],  the  switching 
threshold  in  the  merged  waveguides  will  become  larger  by  the  same  factor  when  expressed  in 
terms  of  power  density  instead  of  y,  because  a  factor  of  8pe  is  embedded  in  the  definition  of  y. 
The  merged  waveguides  will  have  to  be  shorter  by  a  factor  of  8[5e,  because  in  Picture-33,  8  jic 
plays  the  role  of  the  coupling  between  A  i  and  Aj  as  can  be  seen  in  Eq.  (3),  and  because  Lc  must 
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Fig.  5.  The  relative  power  transmission  is  plotted  as  a  function  of  y  for  (a)  Lc  =  L/,  versus 
Lc  =  3 Li,  and  (b)  for  Lc  =  15L/,.  The  NLS  quality  for  Lc  =  3 L/,  and  Lc  =  15L;,  are  lower. 

equal  the  half-beat-length  L/,. 

Injecting  the  double-waveguide  coupler  with  a  complex  beam  profile  of  F\  (x)  to  obtain  the 
A\(Q)/VP  =  1  and  /ti(0)/V/'5  =  0  BCs  may  not  be  possible.  In  many  situations,  it  is  more 
practical  to  inject  the  beam  from  a  waveguide  that  is  identical  to  one  of  the  waveguides  in 
the  multicore  coupler.  For  the  merged  waveguides  (d/ A  =  0),  this  means  that  Fw(x)  is  used 
as  the  injected  and  collected  beam,  instead  of  F\{x).  We  refer  to  this  as  case-2  and  the  NLS 
curve  is  shown  as  dotted  (blue)  in  Fig.  4(b).  Because  the  injected  and  collected  beams  are  not 
properly  mode-matched  to  the  input  and  output  couplers,  respectively,  both  modes  are  excited 
(with  different  powers)  and  we  need  to  use  Aj  (0) /Vp  ss  0.95  and  Ai  (0) / Vp  ss  0.29.  Moreover, 
approximately  1.7%  of  the  injected  power  is  coupled  to  radiation  modes.  The  modulation  depth 
in  the  NLS  curve  in  Fig.  4(b)  for  case-2  in  dotted  (blue)  is  lower  because  of  the  injection  beam 
profile,  resulting  in  a  lower  quality  NLS. 

We  note  that  a  linear  switch  can  fully  transfer  the  power  from  one  waveguide  to  another  if 
Lc  =  nL/,,  where  n  is  an  odd  integer.  Therefore,  all  such  linear  switches  are  equivalent.  However, 
in  the  nonlinear  case  the  quality  of  the  NLS  is  reduced  for  n  >  1  as  shown  in  Fig.  5  and  undesired 
oscillations  are  observed  before  y  is  sufficiently  large  for  full  power  switching. 

3.  Conclusions 

It  is  argued  that  the  switching  power  threshold  in  nonlinear  multimode  junctions  is  larger  than 
in  multicore  junctions,  mainly  because  the  value  of  8fie  is  substantially  larger  in  the  former 
case.  The  XPM/FWM  terms  also  play  a  role  in  setting  the  power  threshold  but  their  impor¬ 
tance  is  orders  of  magnitude  lower  than  that  of  8fie.  When  expressed  in  proper  dimensionless 
parameters,  the  NLS  curves  are  almost  identical  in  a  multimode  and  a  multicore  junction  if  the 
XPM/FWM  terms  are  artificially  switched  off. 

The  presence  of  the  XPM/FWM  terms  do  not  seem  to  play  any  essential  role  in  the  quality 
of  the  nonlinear  switching  curve  of  the  multimode  or  multicore  systems,  especially  in  the  value 
of  the  modulation  depth.  Rather,  injection  of  the  right  combination  of  the  modes  at  the  input 
and  the  collection  of  the  right  combination  of  the  modes  at  the  output  are  the  main  factors 
behind  the  quality  of  the  NLS  curves.  Creating  the  proper  injection  profile  is  often  trivial  in  a 
weakly-coupled  multicore  waveguide;  however,  obtaining  the  right  combination  of  the  modes  in 
a  multimode  junction  is  difficult,  if  not  impossible,  to  achieve.  Therefore,  multimode  junctions 
show  lower  quality  NLS  curves.  This  finding  agrees  with  the  observations  reported  in  Ref.  [8] 
and  provides  an  intuitive  explanation  for  those  results. 

Finite  element  method  has  been  used  to  calculate  the  spatial  profiles  and  propagation  con¬ 
stants  of  the  modes.  Embedded  nonlinear  differential  equation  solvers  in  Mathematic  a  have 
been  used  for  the  NLS  simulations. 
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The  saturable  absorption  characteristics  of  two-,  three-,  and  five-core  one-dimensional  fiber  coupler  arrays  and 
the  seven-core  hexagonal  fiber  coupler  array  are  investigated.  It  is  shown  that  the  performance  of  all  these  satu¬ 
rable  absorbers  is  comparable  and  not  much  is  gained,  if  anything,  by  going  from  a  two-core  nonlinear  coupler 
geometry  to  a  higher  number  of  cores.  This  observation  is  supported  by  the  similarity  of  the  saturable  absorption 
curves,  as  well  as  comparable  pulse  characteristics  obtained  from  the  simulation  of  a  generic  mode-locked  fiber 
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1.  INTRODUCTION 

The  operation  of  a  mode-locked  laser  requires  a  form  of  inten¬ 
sity  discrimination  otherwise  known  as  saturable  absorption 
[1].  Common  saturable  absorbers  (SAs)  include  semiconduc¬ 
tor  SA  mirrors  (SESAMs)  [2],  carbon  nanotubes  [3],  graphene 
[4],  Kerr  lensing  [5],  and  nonlinear  polarization  rotation  (NPR) 
[6].  Desired  attributes  of  SAs  include  fast  response  time 
(relaxation  time),  stability,  long-term  reliability,  ease  of  use, 
appropriate  saturation  fluence,  and  low  loss  [7]. 

For  mode-locked  fiber  laser  applications,  fully  integrable 
SAs  are  desired  in  order  to  take  full  advantage  of  a  robust 
alignment-free  fiber  cavity.  Common  fiber-integrable  SA  solu¬ 
tions  have  deficiencies  that  limit  their  usefulness  for  mode- 
locked  fiber  laser  applications:  SESAMs  generally  have 
picosecond  response  times,  making  them  less  desirable  for 
ultrashort  pulse  generation,  and  their  nonsaturable  losses 
can  result  in  excessive  heating  [8].  The  long-term  reliability 
of  carbon  nanotube-  and  graphene-based  SAs  when  exposed 
to  intense  optical  pulses  is  at  best  questionable  [9],  and  both 
are  limited  in  modulation  depths.  SAs  based  on  NPR  are  very 
sensitive  to  the  slightest  perturbation  of  the  fiber  cavity, 
making  them  unsuitable  outside  the  laboratory. 

Another  problem  with  the  present  SA  technology  concerns 
the  rapid  growth  in  the  maximum  achievable  pulse  energy  and 
peak  power  in  mode-locked  fiber  lasers  over  the  past  decade 
[10-12],  SAs  are  becoming  a  limiting  factor  in  scaling  up  the 
pulse  energy  and  peak  power  to  higher  values.  SESAMs,  car¬ 
bon  nanutubes,  and  graphene  will  likely  not  be  suitable  for 
ultrahigh  pulse  energies  due  to  heat  generation  and  long-term 
reliability.  The  periodic  SA  curve  in  a  NPR  SA  with  respect  to 
the  pulse  power  and  its  sensitivity  to  environmental  perturba¬ 
tions  makes  this  technique  unsuitable  for  ultrahigh-peak- 
power  mode-locked  lasers. 

Nonlinear  mode  coupling  in  waveguide  arrays  (NMCWA), 
including  multicore  fiber  implementations,  offers  an  attractive 
alternative  that  addresses  most  of  the  limitations  faced  by 
common  SAs.  A  dual-core  fiber  laser  geometry  as  an 
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embedded  SA  was  originally  proposed  by  Winful  and  Walton 
[13],  Since  then,  several  theoretical  studies  have  been  con¬ 
ducted  on  nonlinear  mode  coupling  in  semiconductor  coupled 
waveguide  arrays  [14-16],  and  multicore  optical  fibers  [17,18], 
More  recently,  a  seven-core  tapered  optical  fiber  was  ex¬ 
plored  in  the  nonlinear  regime,  and  it  was  shown  to  have 
the  potential  to  be  used  as  an  SA  [19], 

A  mode-locked  fiber  laser  using  an  AlGaAs  waveguide  array 
SA  was  recently  demonstrated  experimentally  [20].  This  dem¬ 
onstration  clearly  serves  as  a  convincing  experimental  proof 
of  concept  for  the  potential  application  of  NMCWA  as  an  SA. 
However,  for  practical  device  applications  and  improved  effi¬ 
ciency,  the  insertion  loss  of  the  waveguide  array  needs  to  im¬ 
prove.  An  all-fiber  SA  in  the  spirit  of  [13,17,18]  seems  like  an 
ideal  solution  to  the  insertion  loss  problem,  because  it  can  be 
easily  spliced  to  other  components  in  the  fiber  cavity. 
Although  semiconductor  arrays  benefit  from  a  larger  nonlin¬ 
ear  coefficient  resulting  in  desirably  short  SAs,  fibers  made 
from  chalcogenide  glasses  with  nonlinear  coefficients  of  up 
to  1000  times  larger  than  silica  can  provide  similar  perfor¬ 
mances  [21,22].  In  addition,  the  nearly  instantaneous  response 
of  nonlinearities  in  optical  fibers  is  ideal  for  ultrashort  pulse 
generation  [22,23]. 

In  this  paper,  we  study  the  performance  of  the  multicore  SAs 
as  a  function  of  the  number  of  waveguide  cores.  Many-core 
geometries  have  been  the  center  of  attention  in  NMCWA 
SAs;  e.g.,  the  NMCWA  studied  in  [14]  has  41  coupled  wave¬ 
guides.  We  explore  and  compare  the  saturable  absorption  char¬ 
acteristics  of  two-,  three-,  and  five-core  one-dimensional  (ID) 
arrays  and  the  seven-core  hexagonal  array  illustrated  in  Fig.  1. 
The  main  question  that  is  addressed  in  this  paper  is  whether 
increasing  the  number  of  cores  improves  the  saturable  absorp¬ 
tion  characteristics  of  these  nonlinear  multicore  arrays. 

The  results  indicate  that  the  performance  of  all  these  SAs  is 
comparable  and  not  much  is  gained,  if  anything,  by  going  from 
a  two-core  nonlinear  coupler  geometry  to  a  higher  number  of 
cores.  This  observation  is  supported  by  the  similarity  of  the 

©  2013  Optical  Society  of  America 


2788  J.  Opt.  Soc.  Am.  B  /  Vol.  30,  No.  11  /  November  2013 


Fig.  1.  Two-,  three-,  and  five-core  ID  arrays  and  the  seven-core  hex¬ 
agonal  array  are  shown  in  a  fiber-optic  geometry.  All  cores  are 
assumed  to  be  identical  in  a  given  geometry.  The  launch  core  is  iden¬ 
tified  with  a  darker  shade. 

saturable  absorption  curves,  as  well  as  comparable  pulse 
characteristics  obtained  from  the  simulation  of  a  generic 
mode-locked  fiber  laser  cavity.  In  the  latter  case,  these  SAs 
are  placed  in  a  fiber  laser  cavity,  in  which  self-starting  stable 
mode-locked  pulses  are  generated.  According  to  the  results, 
one  can  benefit  from  the  simpler  setup  of  a  two-core  fiber  with 
optimized  parameters  to  get  the  desired  output  mode-locked 
pulses  instead  of  using  the  more  complex  multicore  fibers. 

In  the  following,  for  each  multicore  fiber  SA  in  Fig.  1,  the 
optical  beam  is  launched  in  the  launch  core  identified  with  a 
darker  shade  and  is  collected  from  the  same  core  at  the  other 
end  after  propagating  through  the  length  of  the  SA.  The 
dependence  of  the  transmission  through  the  multicore  fiber 
on  the  optical  power  serves  as  the  desired  saturable  absorp¬ 
tion  mechanism. 

2.  DISCUSSION 

In  the  Unear  regime  where  optical  power  is  low,  neighboring 
waveguides  exchange  optical  power  periodically;  the  linear 
coupling  is  caused  by  the  modal  overlap  of  adjacent  wave¬ 
guides  and  is  most  efficient  when  the  adjacent  modes  have 
identical  propagation  constants  [24].  In  Fig.  2,  the  transmis¬ 
sion  through  the  multicore  arrays  of  Fig.  1  is  shown  in  the  lin¬ 
ear  regime  as  a  function  of  the  nonnalized  length  L  relative  to 
c-1;  c  is  the  linear  coupling  coefficient  between  adjacent  cores 
and  is  assumed  to  be  the  same  among  all  neighboring  cores  in 
a  given  geometry.  The  half-beat  length  Lc  is  defined  as  the  dis¬ 
tance  over  which  the  energy  is  maximally  exchanged  from  the 
launch  core  to  the  neighboring  cores  and  depends  on  the  num¬ 
ber  of  fiber  cores  and  the  geometry. 

It  can  be  seen  that  the  half-beat  length  of  a  three-core  fiber 
is  shorter  than  the  two-core  fiber,  because  its  launch  core  is 
coupled  to  “more”  side  cores  with  the  same  strength  and  can 
exchange  the  optical  power  more  efficiently.  The  linear  dy¬ 
namics  in  the  case  of  five  cores  is  a  bit  more  complex  due 


Fig.  2.  Transmission  through  the  two-core  (solid  black),  three-core 
(dashed  red),  and  five-core  (dotted  orange)  fiber  ID  arrays  and  the 
seven-core  hexagonal  array  (dashed-dotted  purple)  is  shown  as  a 
function  of  the  normalized  length  in  the  linear  regime. 
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to  the  differing  power  exchange  efficiencies  of  the  side  cores, 
but  the  first  minimum  happens  for  only  a  slightly  longer  SA 
device.  In  the  case  of  the  seven-core  hexagonal  array,  the  cen¬ 
tral  core  couples  to  six  cores  in  the  outer  ring,  resulting  in  a 
rapid  exchange  of  optical  power  and  shorter  half-beat  length; 
however,  the  power  in  the  central  core  always  remains  larger 
than  zero  where  the  minimum  value  of  transmission  is 
1/7  «  14.3%.  In  all  these  cases,  the  half-beat  length  is  inversely 
proportional  to  the  coupling  strength  and  is  given  by 
cLc  =  £n(n/Z),  where  £n  is  a  numerical  coefficient  for  the 
w-core  scenarios  studied  here.  We  have  £2  =  I,  £:<  = 
1/V2  «  0.71,  £s  =  I/VZI  »  0.77,  and  £7  w  0.378. 

In  the  nonlinear  regime  where  optical  power  is  high, 
self-phase  and  cross-phase  modulation  (SPM  and  XPM)  effects 
alter  the  refractive  index  of  each  waveguide  and  consequently 
detune  the  effective  propagation  constants  of  the  modes, 
reducing  the  power  exchange  efficiency  between  adjacent 
cores  [25].  In  other  words,  if  high  power  is  launched  into  a 
waveguide,  the  effective  couplings  to  neighboring  cores  are  re¬ 
duced  and  the  light  remains  mainly  in  the  launch  core.  If  the 
total  length  of  the  coupler  in  each  case  is  chosen  to  be  equal 
to  the  half-beat  length,  for  an  optical  pulse  transmitting  through 
the  fiber  array,  its  low  intensity  sides  are  efficiently  channeled 
to  the  adjacent  cores,  while  its  high  intensity  center  peak  re¬ 
mains  in  the  launch  core,  resulting  in  a  power-dependent  trans¬ 
mission  and  intensity  discrimination  [13,14],  This  saturable 
absorption  feature  can  be  used  to  produce  the  required  pulse 
shaping  for  stable  and  robust  mode-locked  pulse  trains. 

The  differential  equation  describing  the  propagation  of  the 
electric  field  in  the  nth  core  of  a  multicore  fiber  nonlinear 
coupler,  where  only  the  interactions  between  the  nearest 
neighbors  are  considered,  is  given  by 

^  +  ifn\En\2En  +  ^  X-  CD 

En  is  the  electric  field  envelope  normalized  to  the  pulse  peak 
power  P0  at  the  entrance  of  the  launch  core.  The  n  sums  are 
over  the  cores  adjacent  to  the  nth  core.  The  propagation  dis¬ 
tance  is  rescaled  to  f  =  cz.  yn  =  y*nPa/c  is  a  dimensionless 
parameter  characterizing  the  nonlinear  strength,  where  y% 
is  the  SPM  coefficient  of  the  nth  core  and  /;n  ii  is  the  XPM  co¬ 
efficient  between  the  nth  core  and  the  adjacent  cores.  The 
XPM  terms  are  much  smaller  than  the  SPM  terms  and  are  thus 
neglected,  i.e.,  n„  -n  «  0.  The  dispersion  and  loss  of  the  fiber 
coupler  can  be  neglected  if  the  SA  length  is  sufficiently 
smaller  than  its  dispersion  length  and  effective  length  [23]. 

The  nonlinear  consequences  of  the  coupled  beam  propaga¬ 
tion  in  multicore  arrays  have  been  previously  explored  in  de¬ 
tail,  e.g.,  in  [26],  Here,  we  borrow  from  these  studies  to  obtain 
the  saturable  absorption  curves  of  the  multicore  arrays  of 
Fig.  1.  The  length  of  each  SA  array  is  assumed  to  be  equal 
to  the  half-beat  length  to  ensure  minimum  transmission  in 
the  low-power  linear  regime.  The  saturable  absorption  curves 
are  shown  in  Fig.  3  where  the  transmission  through  the  launch 
core  of  the  multicore  SA  array  is  plotted  as  a  function  of  the 
normalized  nonlinear  strength  y,  which  is  assmned  to  be  the 
same  for  all  cores.  The  saturable  absorption  is  dictated  by  a 
competition  between  the  Unear  coupling  term  versus 

the  nonlinear  SPM  term  iyn\En\2En  in  Eq.  (1),  and  the  transi¬ 
tion  from  low  transmission  to  high  transmission  occurs  when 
these  two  terms  are  comparable. 
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Fig.  3.  Nonlinear  transmission  through  two-core  (solid  black),  three- 
core  (dashed  red),  five-core  (dotted  orange)  fiber  ID  arrays  and  the 
seven-core  hexagonal  array  (dashed-dotted  purple)  is  shown  as  a 
function  of  the  normalized  nonlinear  strength  y,  which  is  assumed 
to  be  the  same  for  all  cores.  The  length  of  each  SA  array  is  chosen 
to  be  equal  to  the  half-beat  length. 


The  three-core  fiber  SA  has  two  adjacent  cores  coupled  to 
its  launch  core  and  consequently  has  a  larger  linear  coupling 
term  compared  with  the  two-core  scenario;  therefore,  more 
optical  power  (larger  SPM  term)  is  required  to  overcome 
the  linear  coupling  term.  That  is  why  the  three-core  SA  re¬ 
quires  slightly  higher  power  for  nonlinear  saturation.  The 
five-core  and  three-core  scenarios  are  nearly  identical,  be¬ 
cause  even  in  the  five-core  scenario,  the  launch  core  is  only 
adjacent  to  two  neighboring  cores.  The  saturable  absorption 
curves  for  a  higher  number  of  cores  in  ID  linear  arrays  are 
practically  identical  to  the  five-core  SA. 

For  the  seven-core  hexagonal  array  SA,  the  nonlinear  SPM 
term  must  overcome  the  linear  coupling  to  six  adjacent  cores 
and  requires  higher  normalized  nonlinearity  to  saturate  the 
SA.  However,  the  SA  curve  closely  resembles  those  of  ID  ar¬ 
rays  and  can  be  mapped  to  the  ID  saturable  absorption  curves 
by  adjusting  the  core-to-core  c  coefficient.  In  fact,  if  SA  devi¬ 
ces  are  chosen  at  the  same  length,  where  c  must  inevitably  be 
smaller  for  the  seven-core  array,  the  saturable  absorption 
curves  for  all  the  geometries  in  Fig.  1  will  be  nearly  identical 
as  a  function  of  the  input  power.  This  can  easily  be  seen  if  the 
“horizontal  scale”  for  the  seven-core  SA  curve  in  Fig.  1  is  re¬ 
scaled  by  a  factor  of  approximately  2/6  to  compensate  for  the 
six  adjacent  cores  versus  the  two  adjacent  cores  in  the  three- 
core  fiber.  The  SA  curves  nearly  overlap  after  such  a  rescaling. 
This  fact  will  result  in  nearly  identical  pulse  characteristics 
in  the  mode-locked  cavity  simulations  of  Fig.  4,  as  will  be 
discussed  shortly. 

There  is  another  notable  difference  in  the  seven-core  SA: 
the  transmission  initially  drops  to  zero  at  finite  y  before  satu¬ 
rating  at  higher  power.  This  effect  is  due  to  the  complex 
coupling  pattern  in  the  seven-core  geometry.  A  nearly  mono¬ 
tonic  increase  can  be  artificially  created  by  shortening  the 
fiber  below  the  half-beat  length  at  the  expense  of  a  smaller 
modulation  depth. 

To  study  the  mode-locking  behavior  of  the  multicore  array, 
it  is  placed  within  the  mode-locking  optical  cavity.  A  typical 
passive  mode-locked  fiber  laser  consists  of  an  active  single¬ 
mode  fiber  with  a  bandwidth-limited  gain,  which  compensates 
for  the  energy  loss  in  the  cavity,  and  an  intensity  discrimina¬ 
tion  device  as  an  SA.  When  the  interactions  of  the  attenuation, 
chromatic  dispersion,  nonlinearity,  and  bandwidth-limited 


L  (mm) 

Fig.  4.  Normalized  peak  power  of  the  mode-locked  pulses  is  shown 
as  a  function  of  the  length  the  SA,  where  the  core-to-core  coupling 
coefficients  c  are  adjusted  so  that  the  length  of  the  SA  is  equal  to 
the  half-beat  length. 


gain  in  the  optical  fiber  cavity  are  combined  with  the  intensity 
discrimination,  a  train  of  mode-locked  pulses  is  formed  from 
white  noise  after  several  numbers  of  round  trips  in  the  cavity. 
These  interactions  are  captured  by  two  sets  of  equations:  one 
is  the  nonlinear  Schrodinger  equation  governing  the  pulse 
propagation  in  the  cavity,  while  the  other  is  Eq.  (1),  describing 
the  saturable  absorption. 

The  pulse  propagation  in  the  cavity  is  given  by  [27] 

fz  =  ^  +  i>m2A  +  9iS)  (*  +  TWjA'  C2) 

where  A(z,  T)  is  the  normalized  electric  field  envelope,  z  is  the 
propagation  distance  along  the  fiber,  and  T  represents  the 
time  measured  in  the  retarded  frame  of  the  pulse  normalized 
by  a  constant  T0,  which  is  the  expected  pulse  duration.  A(z,  T) 
is  normalized  by  the  peak  power  of  the  fundamental  soliton 
Pi s  =  |/?2 1/ where  /?2  is  the  group  velocity  dispersion 
of  the  fiber,  and  y  is  the  fiber  nonlinear  coefficient  [23], 

It  is  also  implicitly  assumed  that  the  wavelength  of  the  pulse 
/l0  =  1.55  pm  falls  within  the  anomalous-dispersion  regime  of 
the  fiber.  The  propagation  distance  z  is  scaled  to  the 
dispersion  length  Ld  =  T\/ \p2\,  and  >/  is  the  dimensionless 
cavity  nonlinear  coefficient  given  by  >/  =  yPlsLD.  The  normal¬ 
ized  attenuation  is  considered  through  a,  while  the  bandwidth- 
limited  gain  of  the  cavity  is  described  by  the  parameters  g(z) 
and  r.  The  dimensionless  parameter  r  controls  the  pulsewidth 
in  the  mode-locking  process  and  is  given  by  r  =  \/(HcoT0)2, 
where  Am  =  (2jic  /  a^)  Ha  is  the  spectral  gain  bandwidth 
and  A/l  is  the  gain  bandwidth  [28],  The  gain  saturation  g(z) 
is  given  by 


g(z)  = 


2ff0 

1  +  \\A\\2/e0  ’ 


(3) 


where  g0  and  e0  are  the  pumping  strength  and  the  cavity 
saturation  energy,  respectively,  and  ||A||2  =  lA^dT  is 
the  pulse  energy. 

The  simulations  are  based  on  Eq.  (2),  while  Eq.  (F)  is  ap¬ 
plied  in  every  round  trip  of  the  laser  cavity.  The  coupling  ef¬ 
ficiency  of  the  cavity  to  the  launch  core  of  the  multicore  SA  is 
assumed  to  be  100%  on  both  the  entrance  and  the  exit  ports; 
therefore,  A  in  Eq.  (2)  becomes  E  at  the  entrance  in  the  launch 
core  of  the  multicore  geometry  in  Eq.  (1)  and  vice  versa  at 
the  exit. 
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For  the  mode-locking  cavity  simulations,  the  initial  pulse 
conditions  are  assumed  to  be  low-amplitude  noise  fluctua¬ 
tions.  The  cavity  length  is  5  m,  while  the  rest  of  the  parameters 
used  are  T0  =  114  fs,  a  =  0.05,  A/1  =  35  nm,  g0  =  0.39,  and 
e0  =  1  [28].  We  also  have  =  -0.0153  ps2/m  resulting  in 
LD  =  0  .84  m.  The  nonlinear  coefficient  of  the  cavity  is  chosen 
as  y  =  1.7  km_1W_1.  The  nonlinear  coefficient  of  the  SA  is 
assumed  to  be  y*  =  0.34  m  1  W1,  which  can  be  obtained 
by  using  a  chalcogenide  glass  fiber  [21], 

Figure  4  shows  the  normalized  power  of  the  mode-locked 
pulses  after  100  round  trips  in  the  laser  cavity,  resulting  in  full 
convergence,  as  a  function  of  the  length  the  SA.  For  each 
length  (and  for  each  SA),  the  core-to-core  coupling  coeffi¬ 
cients  c  are  adjusted  so  that  the  length  of  the  SA  is  equal 
to  the  half-beat  length  for  optimum  performance.  This  choice 
of  the  coupling  coefficients  ensures  that  the  power  is  maxi¬ 
mally  transferred  from  the  launch  core  to  the  other  cores 
in  the  linear  regime.  In  practice,  the  core-to-core  coupling  co¬ 
efficients  c  can  be  designed  and  controlled  via  core-to-core 
separations.  As  we  discussed  above,  the  saturable  absorption 
curves  for  all  multicore  arrays  studied  here  are  nearly  identi¬ 
cal  when  the  SAs  are  of  the  same  length.  Therefore,  it  should 
come  at  no  surprise  that  the  normalized  pulse  peak  power  is 
nearly  the  same  for  all  these  SAs  as  shown  in  Fig.  4. 

In  summary,  the  results  in  Figs.  3  and  4  show  that  nearly 
identical  performance  is  expected  from  the  geometries  that 
are  explored  in  this  manuscript.  While  that  results  are  seem¬ 
ingly  generic  for  popular  SA  geometries,  we  emphasize  that 
this  study  is  not  exhaustive.  An  infinite  number  of  multicore 
geometries  can  be  constructed  even  with  nonuniform  core-to- 
core  couplings  and  nonidentical  core  sizes  applied  to  many 
different  choices  of  laser  parameters.  Therefore,  the  main 
conclusion  one  can  draw  from  this  study  is  that  when  an 
NMCWA  is  to  be  used  as  an  SA,  it  is  important  to  verify 
whether  multiple  cores  yield  any  tangible  performance  bene¬ 
fits  besides  making  the  device  unnecessarily  more  complex. 
Future  efforts  will  focus  on  comparisons  between  multimode 
and  multicore  all-fiber  SA  devices  [29], 
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Abstract:  We  report  on  the  generation  of  new  wavelengths,  mediated  by 
the  stimulated  Raman  scattering  process  and  extending  over  two  octaves 
covering  523  to  1750  nm  wavelength  range,  in  a  standard  telecommunica¬ 
tion  graded-index  multimode  optical  fiber.  Despite  the  highly  multimode 
nature  of  the  pump,  the  Raman  peaks  are  generated  in  specific  modes  of 
the  fiber,  confirming  substantial  beam  cleanup  during  the  stimulated  Raman 
scattering  process. 
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1.  Introduction 


Stimulated  Raman  scattering  (SRS)  is  a  well-known  nonlinear  process  with  numerous  appli¬ 
cations  including  optical  amplifiers,  tunable  lasers,  spectroscopy,  meteorology,  and  optical  co¬ 
herence  tomography.  SRS  was  first  observed  in  silica  glass  fibers  in  1972  by  Stolen  et  al.  [1]. 
They  used  a  frequency-doubled,  pulsed,  Nd:YAG  laser  operating  at  the  532  nm  wavelength  to 
pump  a  single-mode  optical  fiber  and  observed  Stokes  emission  at  545  nm.  Since  then,  several 
groups  have  reported  the  observation  of  SRS  in  optical  fibers  using  various  configurations.  For 
example,  Cohen  and  Lin  [2]  generated  6  cascaded  Raman  peaks  in  a  silica  fiber,  pumped  by  a 
mode-locked,  Q-switched,  Nd:YAG  laser  operating  at  1064  nm.  Rosman  [3]  observed  15  or¬ 
ders  of  cascaded  Raman  peaks  by  pumping  a  silica  fiber  at  532  nm  with  a  frequency-doubled 
Nd:YAG  laser.  Other  configurations  involving  unconventional  fibers  have  been  used  for  ex¬ 
treme  Raman-comb  generation.  For  example,  Couny  et  al.  [4]  demonstrated  the  generation  of 
a  Raman  comb  spanning  wavelengths  from  325  nm  to  2300  nm  in  a  lm-long  hydrogen-filled 
hollow-core  photonic  crystal  fiber. 

In  this  article,  we  report  on  the  generation  of  new  wavelengths,  mediated  by  the  SRS  process 
and  spanning  over  two  octaves,  in  a  standard  graded-index  multimode  fiber  (GIMF).  The  GIMF 
is  pumped  at  a  wavelength  of  523  nm,  and  an  optical-frequency  Raman  comb  is  generated  on  the 
Stokes  side  of  the  pump.  At  high  power  levels  the  measured  spectrum  extends  up  to  1750  nm, 
which  is  the  upper  detection  limit  of  our  optical  spectrum  analyzer  (OSA).  The  generation  of 
such  a  wide  wavelength  range,  covering  two  octaves  from  523  to  1750  nm,  using  a  large-core 
telecommunication-grade  multimode  fiber  distinguishes  our  results  from  those  carried  out  in 
small-core  or  highly  customized  optical  fibers. 

Multimode  optical  fibers  are  easy  to  handle  and  are  also  easy  to  align  to  external  sources; 
however,  their  large  core  diameter  is  perceived  as  undesirable  for  nonlinear  applications.  De¬ 
spite  a  lower  effective  nonlinearity  associated  with  a  larger  core  of  conventional  multimode 
fibers,  the  multimode  nature  of  these  fibers  can  play  an  important  role  in  some  nonlinear  ap¬ 
plications.  In  particular,  the  presence  of  multiple  propagating  modes  with  different  dispersive 
properties  results  in  expanded  phase-matching  opportunities  for  the  generation  of  four-wave 
mixing  (FWM)  signals  in  optical  fibers  [5-7].  The  GIMF  used  in  our  experiments  has  two 
desirable  properties  that  make  it  particularly  suitable  for  SRS  generation.  First,  since  the  ef¬ 
fective  modal  area  of  each  mode  in  the  GIMF  scales  only  as  square  root  of  the  core  area,  the 
effective  nonlinearity  of  some  propagating  modes  is  comparable  with  conventional  single-mode 
fibers  [8].  Second,  a  relatively  high  Ge02  content  in  the  core  of  the  standard  telecommunication 
GIMF  used  in  our  experiments  results  in  a  higher  peak  Raman  gain  coefficient  compared  with 
silica-core  fibers  [9]. 

2.  Generation  of  new  wavelengths  and  spectral  measurements 

The  pump  used  in  our  experiment  is  a  frequency-doubled,  Q-switched,  Nd:  YLF  laser  operating 
at  523  nm  wavelength,  and  its  8-ns  duration  pulses  are  coupled  to  the  input  tip  of  the  fiber 
using  a  microscope  objective.  The  laser  beam  is  not  diffraction-limited,  and  many  modes  in 
the  fiber  are  excited  simultaneously.  The  GIMF  is  a  1-km  long  standard,  50/125-/im,  bare  fiber 
(Corning  ClearCurve  OM2).  The  output  spectrum  is  measured  using  a  CCS200  spectrometer 
(from  Thorlabs)  operating  in  the  range  of  200-1000  nm  and  a  MS9740A  OSA  (from  Anritsu) 
covering  the  600-1750  nm  wavelength  range.  We  measured  the  energy  of  input  pulses  required 
for  reaching  the  Raman  threshold  immediately  after  coupling  into  the  fiber  (less  then  1 -meter 
of  propagation  inside  the  GIMF)  to  be  20.9  fJ.  J;  we  estimate  the  peak  power  to  be  about  2.5  kW. 
The  pulse  energy  decreased  to  about  0.515  /iJ  after  1-km  of  propagation,  which  is  consistent 
with  the  expected  attenuation  of  about  16  dB/km  at  the  pump  wavelength. 

The  sequential  generation  of  cascaded  Raman  peaks  at  four  different  pump  power  levels  is 


Fig.  1 .  Sequential  generation  of  Raman  peaks  at  four  pump  power  levels  as  measured  by  the 
spectrometer.  The  vertical  axis  is  in  arbitrary  “linear”  unit  of  power.  Labels  are  omitted  on 
the  horizontal  wavelength  axis  for  clarity.  The  leftmost  peak  in  each  subplot  is  at  523  nm. 


Fig.  2.  (left)  Cascaded  Raman  peaks  measured  with  the  spectrometer,  (right)  Spectrum 
measured  by  the  OSA  when  pump  power  is  increased  to  just  below  the  burning  threshold 
of  fiber’s  tip;  the  spectral  dip  at  around  1300  nm  and  the  broad  peak  beyond  1400  nm  are 
two  notable  features  in  this  infrared  range. 


shown  in  four  separate  subplots  in  Fig.  1  using  spectrometer  measurements.  The  pump  power 
just  exceeds  the  Raman  threshold  in  subplot  (a)  showing  the  pump  at  523  nm  and  the  first 
low-power  Stokes  line.  As  the  pump  power  is  increased,  the  first  Stokes  line  extracts  power 
from  the  pump  until  it  becomes  strong  enough  to  seed  the  generation  of  next  Stokes  line.  This 
process  continues  and  more  and  more  Raman  peaks  are  generated  with  increasing  pump  power, 
as  shown  in  subplots  (b),  (c),  and  (d)  in  Fig.  1.  The  20  cascaded  peaks  shown  in  Fig.  2  (left) 
extend  from  523  nm  to  just  above  1000  nm  in  wavelength.  The  estimated  input  peak  power  of 
our  pulses  is  22  kW  for  this  figure. 

As  the  pump  power  is  further  increased,  even  more  cascaded  Raman  peaks  appear  beyond  the 
1000  nm  wavelength  range  of  our  spectrometer.  Figure  2  (right)  shows  the  spectrum  measured 
with  the  Anritsu  OSA  in  the  range  of  900-1750  nm  at  the  maximum  pump  power  level  (just 
below  the  burning  threshold  of  fiber’s  input  tip).  The  spectral  dip  at  around  1300  nm  and  the 
broad  peak  beyond  1400  nm  are  two  notable  features  in  this  infrared  range.  In  our  opinion,  the 
broad  peak  centered  at  1600  nm  results  from  the  onset  of  modulation  instability  in  the  presence 
of  anomalous  dispersion.  The  resulting  short  pulses  can  undergo  soliton  fission  and  experience 
Raman-induced  spectral  shifts  as  well  as  soliton-collision-based  spectral  broadening,  resulting 
in  a  broad  supercontinuum  like  feature  [10].  We  stress  that  the  two  plots  in  Fig.  2  should  not  be 
compared  directly  since  they  correspond  to  different  power  levels  and  employ  different  vertical 
scales. 

The  spectral  dip  centered  at  the  1320-nm  wavelength  is  surprising  at  first  sight.  We  think 
its  appearance  is  related  to  a  reduction  in  the  SRS  gain  occurring  near  the  zero-dispersion- 
wavelength  of  the  fiber,  where  the  SRS  gain  is  suppressed  due  to  a  near-perfect  phase-matching 
of  the  FWM  process  [10].  The  dip  at  1320  nm  can  be  seen  more  clearly  in  Fig.  3,  where  we 
plot  the  data  in  Fig.  2  (right)  on  a  logarithmic  power  scale. 
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Fig.  3.  Same  as  Fig.  2  (right)  but  data  plotted  using  a  logarithmic  vertical  scale.  The  input 
pump  power  is  just  below  the  burning  threshold  of  fiber’s  tip. 


Fig.  4.  Raman  comb  measured  with  the  spectrometer  and  plotted  using  frequency  shift  on 
the  horizontal  axis.  Notice  the  presence  of  FWM  peaks  on  the  anti-Stokes  side. 


We  note  that  the  presence  of  efficient  FWM  phase-matching  opportunities  in  the  GIMF  im¬ 
pacts  the  generation  of  the  cascaded  Raman  peaks  [8].  Figure  4  shows  the  spectrometer  data 
in  the  frequency  domain,  by  plotting  on  the  horizontal  axis  the  frequency  shift  of  the  Raman 
comb  relative  to  the  pump  frequency,  and  power  on  a  logarithmic  scale  on  the  vertical  axis. 
Equal  spacing  of  various  peaks  on  the  Stokes  side  is  expected  for  a  cascaded  Raman  process. 
Flowever,  the  presence  of  FWM  frequencies  on  the  anti-Stokes  side  of  the  pump  is  the  most 
notable  feature  in  this  figure.  The  phase-matched  frequency  counterparts  of  these  FWM  idlers 
on  the  Stokes  side  can  affect  the  location  and  amplitude  of  the  cascaded  Raman  peaks.  We 
observed  that  the  strength  of  the  FWM  signal  depended  on  launch  conditions,  and  FWM  was 
absent  (or  highly  suppressed)  in  some  of  our  measurements. 

In  order  to  explore  this  effect,  we  slightly  offset  the  input  pump  beam  so  that  beam  center  is 
not  aligned  with  the  center  of  the  fiber  core.  The  pump  now  excites  the  GIMF  modes  with  dif¬ 
ferent  power  ratios,  resulting  in  efficient  FWM  in  a  different  set  of  phase-matched  wavelengths. 
The  result  is  a  shift  in  the  position  of  the  spectral  combs.  The  red  and  green  spectra  in  Fig.  5  are 
measured  before  and  after  offsetting  the  pump  laser,  respectively.  The  shift  seems  to  be  seeded 
around  the  third  cascaded  Raman  peak,  separated  by  about  50  THz  from  the  pump  frequency, 
which  is  also  consistent  with  the  location  of  a  FWM  peak  in  Fig.  4.  Similar  observations  of  the 
effect  of  the  FWM  processes  on  shifting  the  SRS  spectrum  have  been  reported  by  Sharma  et 
al.  [11];  they  have  shown  that  the  cascaded  Raman  peaks  can  shift  depending  on  which  modes 
are  excited  by  the  pump  laser. 

3.  Spatial  beam  profiles  and  beam  cleanup  in  SRS  generation 

We  also  measured  the  transverse  intensity  profile  of  the  output  beam  by  a  CCD  camera  using 
several  bandpass  color  filters  with  a  bandwidth  of  about  10  nm  (full  width  at  half-maximum). 
The  results  are  shown  in  Fig.  6.  The  leftmost  profile  is  measured  without  a  color  filter,  and 
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Fig.  5.  SRS  spectra  measured  before  (red  solid  line)  and  after  (dashed  green  line)  a  slight 
offset  of  the  pump  beam  from  the  fiber  center.  Shift  of  Raman  peaks  is  caused  by  different 
FWM  conditions  caused  by  excitation  of  different  fiber  modes. 


Fig.  6.  Measured  spatial  profiles  using  a  CCD-based  beam  profiler.  The  leftmost  image  is 
measured  with  no  filters.  The  other  4  images  are  obtained  using  color  filters  centered  at 
610  nm  (second  from  left),  700  nm  (third  from  left),  770  nm  (fourth  from  left),  and  890  nm 
(rightmost). 


the  interference  of  multiple  modes  can  be  clearly  observed  as  a  speckle  pattern.  When  a  color 
filter  centered  at  610  nm  is  placed  in  front  of  the  beam  (second  from  left),  we  observe  a  narrow 
round  spot  that  appears  to  correspond  to  the  spatial  profile  of  the  fundamental  LGoo  mode  of 
our  GIMF  (LG  stands  for  Laguerre -Gauss).  A  donut-shape  spot  (third  from  left)  is  observed 
when  a  color  filter  centered  at  700  nm  is  placed  in  front  of  the  beam;  it  may  correspond  to  the 
LGoi  mode  of  our  fiber. 

Spatial  beam  profiles  corresponding  to  higher-order  fiber  modes  were  also  seen  in  our  ex¬ 
periment.  As  two  examples,  the  two  rightmost  images  of  Fig.  6  show  profiles  corresponding 
to  the  LGio  and  LG20  modes.  They  were  obtained  by  using  optical  filters  with  a  10-nm  pass- 
band  centered  at  770  nm  and  890  nm  wavelengths,  respectively.  The  most  notable  feature  we 
want  to  stress  is  that  a  GIMF  can  be  used  as  a  device  that  not  only  shifts  the  pump  wavelength 
toward  the  red  side  through  SRS  but  also  performs  the  beam  cleanup  owing  to  the  fact  that 
different-order  Raman  peaks  generally  propagate  in  different  modes  of  the  fiber.  An  explana¬ 
tion  of  SRS-induced  beam  cleanup  in  graded-index  multimode  optical  fibers  can  be  found  in 
Ref.  [12].  Chiang  reported  similar  results  for  higher  order  SRS  combs  in  a  30-m-long  fiber  [13]. 
However,  only  LP01  mode  (corresponding  to  LGoo  mode  here)  was  observed  for  a  1-km-long 
fiber.  In  our  experiments,  we  observed  higher-order  modes  even  for  a  1-km-long  fiber,  and  beam 
cleanup  was  not  at  the  same  level  reported  in  Ref.  [13]. 

4.  Conclusion 

We  have  used  a  standard,  telecommunication-grade,  graded-index  multimode  fiber  for  SRS 
generation  by  pumping  it  at  523  nm  with  8-ns  pulses.  We  observe  multiple  cascaded  Raman 
peaks  extending  up  to  1300  nm.  Beyond  that  wavelength,  the  nature  of  dispersion  changes 
from  normal  to  anomalous  since  our  fiber  has  its  zero-dispersion  wavelength  near  1320  nm.  At 


higher  pump  powers,  in  addition  to  the  multiple  cascaded  Raman  peaks,  we  observe  a  single 
broadband  spectral  peak,  extending  from  1350  to  1750  nm.  Its  origin  appears  to  be  related  to 
the  formation  of  solitons  through  modulation  instability  and  their  subsequent  fission  induced 
by  higher-order  dispersion  and  intrapulse  raman  scattering.  Such  features  have  been  observed 
in  the  past  for  single-mode  fibers.  Our  experiments  show  that  a  supercontinuum  can  also  form 
in  a  multimode  fiber.  The  multimode  nature  of  the  fiber  can  also  be  useful  from  a  practical 
standpoint.  For  example,  we  observed  that  different  spectral  peaks  have  spatial  patterns  that 
correspond  to  different  fiber  modes.  This  feature  can  be  useful  for  beam  cleanup.  Future  efforts 
will  focus  on  extending  the  spectrum  to  the  infrared  region  and  on  stabilizing  the  frequency  and 
power  of  individual  comb  lines  for  practical  applications. 
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We  propose  an  all-fiber  mode-locking  device,  which  operates  based  on  nonlinear  switching  in  a  novel  two- 
concentric-core  fiber  structure.  The  design  is  particularly  attractive  given  the  ease  of  fabrication  and  coupling 
to  other  components  in  a  mode-locked  fiber  laser  cavity.  The  nonlinear  switching  in  this  coupler  is  studied, 
and  the  relative  power  transmission  is  obtained.  The  analysis  shows  that  this  nonlinear  switch  is  practical  for 
mode-locking  fiber  lasers  and  is  forgiving  to  fabrication  errors.  ©  2014  Optical  Society  of  America 
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A  nonlinear  optical  device,  referred  to  as  a  saturable 
absorber  (SA),  is  required  to  start  the  lasing  process  in 
passively  mode-locked  lasers.  The  power-dependent 
transmission  of  SAs  creates  an  intensity  contrast  in  the 
laser  cavity,  which  leads  to  the  suppression  of  random 
fluctuations  and  the  generation  of  optical  pulse  trains. 
Common  fiber-integrable  SAs  have  deficiencies  such 
as  long-term  reliability,  thermal  issues,  and  sensitivity 
to  environmental  perturbations,  which  limit  their  useful¬ 
ness  for  mode-locked  fiber  laser  applications.  Also, 
common  SAs  have  become  a  limiting  factor  in  boosting 
the  pulse  energy  and  peak  power  of  mode-locked  lasers 
to  higher  values.  Hence,  there  is  a  strong  need  to  develop 
robust  alternative  mode-locking  techniques. 

Recently,  nonlinear  mode  coupling  in  waveguide  ar¬ 
rays  and  multicore  fibers  [1-5]  and  also  nonlinear  multi¬ 
modal  interference  in  graded-index  multimode  fibers  [6] 
have  been  introduced  as  alternative  techniques  that 
address  most  of  the  limitations  faced  by  common  SAs. 
In  Ref.  [7],  we  compared  the  nonlinear  switching  behav¬ 
ior  in  multicore  versus  multimode  fibers  and  showed  that 
a  much  higher  switching  power  is  required  in  the  latter 
case,  due  to  the  large  difference  in  the  propagation 
constant  of  the  modes.  This  feature  makes  multimode 
fiber  couplers  potential  mode-locking  devices  for  future 
high-peak-power  lasers.  For  lower  peak  powers,  as  is 
the  case  in  most  present-day  mode-locked  fiber  lasers, 
it  is  preferred  to  use  multicore  couplers. 

Our  recent  study  on  multicore  fiber  array  couplers 
showed  that  going  from  a  two-core  fiber  geometry  to  a 
higher  number  of  cores  does  not  improve  the  nonlinear 
switching  performance  of  the  device  considerably,  if  any 
at  all  [8].  Hence,  according  to  our  studies,  a  two-core  fi¬ 
ber  coupler  with  optimized  parameters  is  the  optimum 
solution  to  achieve  an  all-fiber  mode-locking  device. 

In  this  Letter,  we  propose  a  novel  all-fiber  nonlinear 
switching  and  mode-locking  device  in  the  form  of  a  tri¬ 
ple-clad  fiber,  as  shown  in  Fig.  3.,  where  the  refractive 
index  of  the  core  and  that  of  the  second  cladding  layer 
are  the  same.  We  will  refer  to  this  structure  as  a  two- 
concentric-core  (TCC)  fiber,  which  is  composed  of  a 
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central  circular  core  and  a  surrounding  ring  core 
[9,10].  The  main  reason  for  choosing  this  geometry  is 
its  simple  fabrication  process.  Given  that  conventional 
fibers  are  fabricated  radially,  maintaining  the  sensitive 
specifications  of  this  fiber  during  the  drawing  process 
is  easier  than  with  nonconcentric  conventional  two-core 
fibers  [IT].  Another  related  reason  is  that  aligning  the 
central  core  of  this  fiber  with  the  single-mode  fibers 
(SMFs)  in  the  laser  cavity  is  much  easier  and  more  robust 
compared  with  conventional  two-core  geometries. 

In  this  Letter,  the  parameters  of  this  two-core-fiber 
mode-locking  device  are  designed,  and  the  nonlinear 
propagation  of  a  pulse  in  this  coupler  is  analyzed  numeri¬ 
cally.  To  lower  the  required  switching  power,  it  would  be 
preferred  to  build  the  coupler  from  a  material  with  a 
higher  nonlinearity  compared  to  silica  [12].  Chalcogenide 
glass  is  a  suitable  candidate,  because  it  has  a  nonlinear 
coefficient  up  to  1000  times  larger  than  silica  [13], 

The  all-fiber  mode-locking  device  is  shown  in  Fig.  2, 
where  the  beam  is  injected  into  the  central  core  of  the 
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Fig.  1.  Cross-sectional  view  and  the  refractive  index  profile  of 
the  fiber. 
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Fig.  2.  All-fiber  nonlinear  switching  device  composed  of  a 
TCC  coupler  placed  in  between  two  SMFs. 


TCC  fiber  coupler  using  an  input  SMF  and  is  collected  at 
the  output  from  another  SMF.  The  injected  beam  from 
the  input  SMF  is  in  the  form  of  an  azimuthal  symmetric 
nearly  Gaussian  beam,  so  by  centrally  aligning  the  core  of 
the  input  SMF  with  the  central  core  of  the  TCC  coupler, 
only  zero-angular-momentum  modes  will  be  excited  in 
the  coupler  [14].  For  low-input  optical  powers  launched 
into  the  nonlinear  coupler,  optical  power  is  periodically 
exchanged  among  the  two  cores  because  of  their  modal 
overlap.  To  achieve  maximum  linear  coupling,  the 
parameters  of  the  fiber  coupler  should  be  designed  such 
that  the  two  cores  are  phase  matched  [15].  For  high-input 
optical  powers,  nonlinear  effects  arise  and  alter  the  re¬ 
fractive  index  of  each  waveguide.  This  will  detune  the 
effective  propagation  constants  of  the  modes  and  will 
consequently  reduce  the  power  exchange  efficiency  be¬ 
tween  the  cores,  which  means  that  light  will  be  mainly 
preserved  in  the  launch  core  [16],  Hence,  if  an  optical 
pulse  is  transmitted  through  a  fiber  coupler  with  a  length 
equal  to  the  coupler’s  half-beat  length,  its  low-intensity 
tails  will  be  linearly  coupled  to  the  adjacent  core,  while 
its  high-intensity  center  peak  will  remain  in  the  launch 
core  [17],  If  only  the  light  in  the  launch  core  is  collected 
at  the  output,  the  required  power-dependent  transmis¬ 
sion  for  mode-locking  is  obtained  [1,2], 

The  fiber  profile  is  designed  such  that  it  supports  only 
two  zero-angular-momentum  guided  modes,  referred  to 
as  the  first-  and  second-order  guided  modes.  The  electric 
field  in  the  TCC  can  be  expressed  as  the  summation  of 
these  guided  modes  as 


E{x ,  y ,  z,  t)  =  (z,  t)Ffl(x ,  y)el(-03 at  ,  (1) 

p= 1,2 


where  A^,(z,  t)  is  the  slowly  varying  envelope  of  the 
electric  field  of  the  /,/th  mode  with  the  normalized  spatial 
distribution  F^(x^y)  and  propagation  constant  /10/,.  The 
generalized  nonlinear  Schrodinger  equation  (GNLSE) 
describing  the  longitudinal  evolution  of  the  pulse  in  this 
fiber  coupler  can  be  written  as  [18,19] 
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The  indices  can  take  the  value  of  1  or  2,  corresponding 
to  first-  and  second-order  modes,  a  is  the  fiber  attenua¬ 
tion,  h  is  the  group  delay  per  unit  length,  and  p2,P  is 
the  group  velocity  dispersion  (GVD)  parameter  of  the 
//th  mode.  We  define  =  -d/I02  =  /?0,i  -  A),ref  and 
<%. i  =  -Sfl  1,2  =  Pi,i  -/?i,ref>  where  =  (fin  +  /?i,2)/2 
for  i  =  0, 1.  n2  is  the  nonlinear  index  coefficient,  and 


w(}  is  the  carrier  frequency.  The  nonlinear  coupling  coef¬ 
ficients  are  given  by 

fpvKp  =  /  drd  yF;FvFKF*p,  (3) 

where  the  mutually  orthogonal  spatial  profiles  are  as¬ 
sumed  to  be  normalized  according  to  ff  F2dxdy  = 
ff  Fidzdy  =  1.  The  spatial  mode  profiles  of  the  guided 
modes  of  the  TCC  coupler,  Fl  and  F2,  and  also  that  of 
the  center  core,  Fc,  can  be  found  using  the  finite  element 
method  (FEM)  [20],  where  Fc  is  calculated  for  the  central 
core  in  the  absence  of  the  outer  ring.  The  beam  is 
launched  into  and  collected  from  the  center  core  of  the 
TCC  fiber  coupler.  The  transmittance  is  defined  as  the 
energy  at  the  output  of  the  center  core  divided  by  the  to¬ 
tal  energy  injected  into  the  coupler,  given  by 

_  =  /-+q0°°  \AC(LC,  t)\2dt 

/-+“(l^i((U)|2  +  \a2(0,  t)\2)dt’ 

where  Lc  is  the  half-beat  length,  Pout  =  |AC|2  is  the  col¬ 
lected  output  power  from  the  coupler,  and  Ac  is  the  field 
amplitude  of  the  center  core,  defined  as 


Ac(z,t)  =Al(z,t )  //  F*F1dxdy  +  A2(z,t)  //  F*F2dxdy 


(5) 


In  the  following  simulations,  we  have  assumed  that  the 
fiber  material  is  As2Se3  chalcogenide  glass.  The  operat¬ 
ing  wavelength  is  considered  a{)  =  1550  ran.  The  refrac¬ 
tive  index  of  both  cores  is  nco  =  2.834,  while  that  of  the 
cladding  is  nc\  =  2.831.  The  central  core  has  a  radius  of 
r1  =  4.1  pm,  and  the  inner  radius  of  the  ring  is 
r2  =  8.2  pm,  while  its  width  is  3.94  pm.  The  nonlinear  in¬ 
dex  is  n2  =  1 . 1  x  10-17  m2/W,  and  the  fiber  propagation 
loss  is  a  =  1  dB/rn  1 13,21].  The  dispersion  of  the  funda¬ 
mental  mode  is  equal  to  D  =  -826  ps/nm/km,  while  that 
of  the  second-order  mode  is  D  =  -819  ps/nm/km.  These 
values  include  both  the  material  and  the  waveguide 
dispersion.  The  material  dispersion  is  calculated  using 
the  Sellmeier  equation  [22]  and  then  incorporated  in 
the  FEM  to  evaluate  the  total  chromatic  dispersion. 
The  fiber  exhibits  normal  dispersion  at  1550  nm.  It  should 
be  noted  that  fibers  with  such  a  layered  refractive  index 
profile  are  frequently  used  for  dispersion  compensation 
due  to  their  large  negative  dispersion  values  [23]. 

With  these  parameters,  each  individual  core  supports 
only  1  zero-angular-momentum  mode.  It  should  be 
pointed  out  that  in  order  to  minimize  the  splicing  losses 
between  the  silica  SMF-28  fibers  of  the  cavity  and  the 
TCC,  the  dimensions  and  numerical  aperture  of  the  cen¬ 
tral  core  are  chosen  to  be  similar  to  those  of  the  cavity 
SMFs.  The  Fresnel  reflection  loss,  due  to  the  refractive 
index  difference  between  silica  and  chalcogenide  glass, 
is  around  0.45  dB  on  each  end  of  the  TCC  fiber,  which 
can  be  potentially  suppressed  by  using  antireflection 
coating  [24].  Whether  the  facet  reflection  can  be  toler¬ 
ated  or  would  require  impedance  matching  depends  on 
the  configuration  of  the  mode-locked  fiber  laser  cavity, 
and  is  beyond  the  scope  of  this  analysis.  The  main 
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advantage  of  this  to  silica-based  structures,  as  will  be 
shown  later,  is  its  substantially  lower  switching  power. 
In  what  follows,  we  will  refer  to  the  power  after  reflec¬ 
tion  and  coupling  losses  at  the  input  end  as  the  input 
power,  and  the  power  before  reflection  and  coupling 
losses  at  the  output  end  as  the  output  power.  The  length 
of  the  TCC  coupler  is  1.74  mm,  given  by  Lc  =  n/(2C), 
where  C  =  d/i0]  is  the  linear  coupling  coefficient  among 
the  two  propagating  modes  [151. 

We  assume  that  a  hyperbolic  secant  optical  pulse  in  the 
form  of /?(().  t)  =  P\/z  sech(£/£0)  is  launched  into  the  non¬ 
linear  coupler  through  the  input  SMF,  where  P0  is  the 
peak  power,  and  the  full  width  at  half-maximum  (FWHM) 
pulse  duration  is  T 0  »  1.763£0.  The  nearly  Gaussian  spa¬ 
tial  mode  profile  of  the  input  beam  with  a  temporal 
power  profile  of  Pc(t )  excites  the  guided  modes  of  the 
TCC  fiber  with  power  profiles  P\(t)  and  P-At)  with  nor¬ 
malized  peak  powers  p1  and p2,  respectively,  where  it  can 
be  shown  that  pi  +  p2  ~  1-  The  normalized  peak  excita¬ 
tion  power  of  each  guided  mode  at  the  TCC  input  facet  is 
calculated  using  the  overlap  integral  of  the  input  beam 
and  that  particular  mode,  where  in  this  case  it  is  found 
that  Pi  =  62.6%  and  p2  =  37.2%. 

We  have  applied  a  split-step  Fourier  method  [251  to 
solve  Eq.  (2)  and  analyze  the  propagation  of  the  pulse 
inside  the  nonlinear  coupler.  For  low  powers,  the  TCC 
fiber  operates  similarly  to  a  linear  directional  coupler 
at  half-beat  length,  and  most  of  the  pulse  energy  is  trans¬ 
ferred  to  the  ring  core.  For  intense  pulses,  the  central 
portion  of  the  pulse  remains  in  the  central  waveguide, 
while  the  low-intensity  wings  of  the  pulse  are  transferred 
to  the  ring  core.  In  Fig.  3,  the  normalized  input  pulse 
power  launched  into  the  central  core  and  that  of  the  out¬ 
put  pulse  are  shown  for  a  FWHM  pulse  duration  of 
T0  =  1  ns.  The  reason  for  choosing  a  nanosecond  pulse 
width  is  that  in  most  modern  lasers  the  SA  is  a  critical 
device  only  for  allowing  the  pulse  to  build  up  from  the 
noise  and  start  the  lasing  process,  while  a  narrow  spec¬ 
tral  filter  is  responsible  for  the  pulse  shaping  [26].  Hence, 
nanosecond  pulses  which  are  pertinent  to  the  starting 
process  are  considered  here. 

The  effect  of  the  pulse  duration  on  the  switching  dy¬ 
namics  can  be  seen  in  Fig.  4,  where  we  have  calculated 
the  relative  power  transmission  as  a  function  of  the  input 


Fig.  3.  Normalized  input  and  output  pulse  shapes  for  a  pulse 
width  of  T0  =  1.763£0  =  1  ns.  The  input  peak  power  is  5  kW. 


peak  power  for  various  input  pulse  durations.  The  rela¬ 
tive  power  transmission  for  pulse  widths  larger  than  1  ps 
are  the  same  as  that  of  T0  =  1  ps.  It  can  be  observed  that 
a  higher  input  power  is  required  to  switch  the  output 
from  rniin  to  rmax  as  the  duration  of  the  pulse  gets  shorter. 
Also,  the  maximum  achievable  transmission  decreases  as 
the  pulse  duration  is  shortened.  The  reason  is  that  the 
dispersion  effects  become  more  significant  for  shorter 
pulse  widths,  and  the  resultant  pulse  broadening  pre¬ 
vents  complete  power  transfer.  It  should  be  noted  that, 
as  the  input  peak  power  is  increased,  nonlinear  terms  be¬ 
come  larger,  and  because  of  the  large  amount  of  the 
induced  frequency  chirp  imposed  on  the  pulse  due  to 
self-phase  modulation,  even  weak  dispersive  effects  lead 
to  significant  pulse  shaping.  Hence,  for  high  input 
powers,  the  output  pulse  would  become  rectangularly 
shaped  with  relatively  sharp  leading  and  trailing  edges 
[25].  This  limits  the  upper  level  of  the  input  peak  power 
to  these  couplers.  The  damage  threshold  of  chalcogenide 
glasses  should  also  be  considered  while  designing  these 
couplers.  The  optical  damage  mechanisms  for  these  fi¬ 
bers  are  not  well  known;  however,  the  obtained  switch¬ 
ing  power  in  our  design  is  not  far  from  the  damage 
threshold  of  another  chalcogenide  glass  composition, 
As2S3  [27,28]. 

The  switching  power  is  considered  to  be  the  power  at 
which  the  output  from  the  central  core  contains  50%  of 
the  input  energy  injected  into  the  fiber.  While  this  power 
is  slightly  above  3  kW  for  the  chalcogenide  fiber  coupler 
considered  here,  it  is  of  the  order  of  MWs  for  silica  fibers 
with  the  same  structure  [29].  Also,  one  can  reduce  the 
switching  power  by  increasing  the  distance  between 
the  two  cores  and  reducing  the  coupling  between  them. 
For  example,  increasing  the  gap  between  the  two  cores 
to  9  pm  results  in  a  switching  power  of  800  W.  However, 
the  length  of  the  coupler  is  longer  in  this  case,  and  the 
strong  dispersive  effects  overcome  the  nonlinear  pulse 
shortening  mechanism  if  the  pulse  is  too  short. 

During  the  fabrication  process,  random  core  diameter 
fluctuations  will  be  inevitable;  the  effect  of  these  fabrica¬ 
tion  imperfections  on  the  performance  of  the  coupler  is 
studied  in  Ref.  [29],  and  it  is  shown  that  these  imperfec¬ 
tions  are  tolerable  as  long  as  the  changes  they  impose  on 
the  propagation  constant  of  the  modes  are  smaller  than 
the  coupling  between  the  two  cores.  For  the  parameters 


Input  peak  power  (kW) 

Fig.  4.  Relative  power  transmission  is  plotted  as  a  function  of 
the  input  peak  power  for  different  pulse  widths. 
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considered  here,  typical  ±2%  core  size  variations  are 
tolerable. 

In  summary,  we  have  proposed  a  novel  concentric- 
core  fiber  structure  made  of  chalcogenide  glass  to  be 
used  as  an  all-fiber  mode-locking  element.  The  nonlinear 
switching  behavior  of  this  coupler  is  analyzed  by  carrying 
out  a  normal  mode  analysis  of  the  waveguide.  The  chal¬ 
cogenide  glass  has  a  large  chromatic  dispersion.  In  order 
to  avoid  substantial  pulse  broadening  due  to  dispersion, 
the  length  of  the  coupler  (half-beat  length)  should  be  re¬ 
duced  as  much  as  possible.  This,  in  turn,  can  be  achieved 
by  reducing  the  separation  between  the  central  core  and 
the  ring  core  and  increasing  the  coupling  between  them. 
On  the  other  hand,  increasing  the  coupling  between  the 
cores  increases  the  switching  power,  which  may  not  be 
desirable.  Therefore,  the  optimal  design  is  a  compromise 
between  the  pulse  broadening  and  switching  power.  The 
pulse  broadening  and  switching  power  can  also  be  modi¬ 
fied  by  changing  the  geometrical  parameters  and  glass 
compositions. 
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We  explore  the  practical  challenges  that  should  be  addressed  when  designing  a  multicore  fiber  coupler  for  non¬ 
linear  switching  or  mode-locking  applications.  The  inevitable  geometric  imperfections  formed  in  these  fiber  cou¬ 
plers  during  the  fabrication  process  affect  the  performance  characteristics  of  the  nonlinear  switching  device. 
Fabrication  uncertainties  are  tolerable  as  long  as  the  changes  they  impose  on  the  propagation  constant  of  the 
modes  are  smaller  than  the  linear  coupling  between  the  cores.  It  is  possible  to  reduce  the  effect  of  the  propagation 
constant  variations  by  bringing  the  cores  closer  to  each  other,  thereby  increasing  the  coupling.  However,  higher 
coupling  translates  into  a  higher  switching  power,  which  may  not  be  desirable  in  some  practical  situations.  There¬ 
fore,  fabrication  errors  limit  the  minimum  achievable  switching  power  in  nonlinear  couplers.  ©  2014  Optical 
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1.  INTRODUCTION 

Nonlinear  mode  coupling  in  multicore  fiber  couplers  has  re¬ 
cently  attracted  considerable  attention  owing  to  their  applica¬ 
tions  in  nonlinear  switching,  and  in  particular,  in  all-optical 
mode  locking  of  fiber  lasers  [1—6] .  A  multicore  fiber  operates 
as  a  linear  directional  coupler  at  low  optical  powers,  where 
power  is  exchanged  periodically  among  neighboring  wave¬ 
guides.  Linear  coupling  is  most  efficient  when  the  cores  are 
phase  matched,  or  in  other  words,  their  propagating  modes 
have  identical  propagation  constants  [7].  At  higher  optical 
powers,  the  nonlinear  refractive  indices  of  the  cores  are 
altered  due  to  nonlinear  effects.  Consequently,  the  effective 
propagation  constants  of  the  participating  modes  are  detuned, 
which  reduces  the  power  exchange  efficiency  between  the 
neighboring  cores  and  as  a  result,  the  majority  of  light  remains 
in  the  launch  core  [8].  Therefore,  the  different  behaviors  that  a 
multicore  fiber  shows  based  on  the  intensity  of  the  input 
optical  power  make  it  a  suitable  device  for  nonlinear  switch¬ 
ing  and  mode-locking  applications. 

Nonlinear  switching  and  mode-locking  schemes,  which  rely 
on  nonlinear  coupling  among  the  propagating  modes,  deal 
with  practical  constraints.  For  instance,  the  need  for  a  high 
nonlinear  index  of  refraction  limits  the  application  of  long- 
period  fiber  gratings  in  nonlinear  switching  [9].  In  dual-core 
fiber  ring  lasers,  preserving  strong  resonant  coupling  between 
the  core  modes  over  the  required  cavity  lengths  is  difficult 
[10];  similarly,  maintaining  the  coupling  among  the  modes 
of  a  graded-index  multimode  fiber  over  long  lengths  is  a  chal¬ 
lenge  [11].  In  multicore  fiber  couplers,  uncertainties  such  as 
fluctuations  in  the  core  radii  happen  both  in  the  preform  fab¬ 
rication  and  the  drawing  of  the  fiber.  These  errors  affect  the 
core-to-core  couplings.  In  some  applications,  such  as  image 
transport,  errors  and  variations  are  incorporated  intentionally 
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in  the  cores  to  enhance  the  image  contrast,  because  efficient 
core-to-core  coupling  can  be  detrimental  to  the  image  trans¬ 
port  quality  [12,13],  By  contrast,  applications  such  as  nonlin¬ 
ear  switching  require  accurate  core-to-core  coupling,  so, 
fabrication  errors  degrade  the  performance  of  multicore 
nonlinear  fiber  couplers. 

In  this  paper,  to  get  a  realistic  understanding  of  the  perfor¬ 
mance  of  multicore  fiber  couplers  in  the  presence  of  fabrica¬ 
tion  uncertainties,  the  limitations  that  these  errors  impose  on 
the  performance  of  a  nonlinear  switch  are  investigated  at 
three  levels:  (1)  effects  on  the  propagation  constant  of  the 
propagating  modes  [1,  relative  to  the  core-to-core  coupling; 
(2)  effects  on  the  first  derivative  of  f),  otherwise  known  as 
the  group  velocity;  and  (3)  effects  on  the  second  derivative 
of  [),  known  as  the  group  velocity  dispersion  (GVD).  These 
issues  will  be  addressed  in  the  following  sections  of  the  paper. 

2.  UNCERTAINTIES  IN  FIBER 
FABRICATION 

In  the  linear  regime,  light  propagation  in  each  core  of  a  two- 
core  fiber  coupler  is  generally  described  using  the  coupled 
mode  theory,  and  can  be  written  as 

“■(!;)=(<?:  A)(t)-BA-  « 

where  Aa  and  Ah  are  the  amplitudes  of  the  modes  in  each  indi¬ 
vidual  waveguide.  It  is  assumed  that  the  propagating  modes  in 
both  waveguides  have  identical  propagation  constants,  fi0, 
and  coupling  coefficients,  Cab  =  Cba.  Instead  of  considering 
the  modes  of  individual  waveguides  and  their  couplings, 
an  alternative  approach  is  to  study  the  propagation  of  the 
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fundamental  modes  of  the  complete  index  profile  of  the  cou¬ 
pler,  known  as  the  symmetric  and  antisymmetric  supermodes. 
The  spatial  mode  profile  of  the  supermodes,  Fh(x,y),  can  be 
approximated  as  a  superposition  of  the  guided  modes  of  each 
individual  core  as  [14,15] 

Fh(x,y)Ki  22  y ) ’  (2) 

/i=  1,2 

where  Ff,  (x.  y)  is  the  mode  profile  of  the  fundamental  guided 
mode  propagating  in  each  core  and  kt,  is  the  weight  factor  of 
the  pth  core.  The  eigenvectors  and  eigenvalues  of  B  provide 
the  weight  factors  and  propagation  constants  of  the  superm¬ 
odes,  respectively.  The  eigenvectors  of  B,  which  are  (})  and 
l-,1),  lead  to  equal  excitation  of  both  supermodes  at  the  fiber 
input.  However,  if  due  to  fabrication  uncertainties  the  propa¬ 
gation  constant  of  one  waveguide  changes  slightly  to  p0  +  dfl, 
the  eigenvectors  of  B  change  as  well,  and  consequently,  the 
input  beam  is  coupled  to  the  symmetric  and  antisymmetric 
supermodes  unequally.  This  reduces  the  power  exchange 
efficiency  between  the  two  cores  and  degrades  the  nonlinear 
switching  performance  of  the  coupler  by  lowering  the  modu¬ 
lation  depth  of  the  mode-locking  device,  as  also  seen  in 
Ref.  [15].  In  order  to  avoid  this  and  observe  perfect  linear 
switching  for  low  powers,  the  fluctuations  in  [St]  should  be 
much  smaller  than  the  coupling  between  the  propagating 
modes.  Hence,  we  can  define  8p  «  Cab  as  the  cutoff  for  the 
acceptable  range  of  fabrication  errors;  in  this  case,  around 
70%  of  the  power  is  in  one  supermode,  whereas  the  other 
one  carries  only  30%. 

To  present  the  aforementioned  arguments,  we  explore  the 
nonlinear  switching  behavior  of  a  two-concentric-core  (TCC) 
fiber  geometry.  This  fiber  coupler  has  TCCs,  consisting  of  a 
central  circular  core  and  a  surrounding  ring  core,  as  shown 
in  Fig.  1.  The  refractive  index  profile  of  the  fiber  considered 
can  also  be  seen  in  this  figure.  Owing  to  its  cylindrical 


Fig.  1.  Cross-sectional  view  and  refractive  index  profile  of  the  fiber 
structure. 


symmetric  configuration,  this  fiber  can  be  fabricated  using 
the  same  techniques  used  for  conventional  fibers,  which  is 
more  straightforward  compared  with  nonconcentric  two-core 
fibers.  This  simple  fabrication  process  is  the  main  reason  for 
choosing  this  specific  geometry.  One  other  reason  is  that  the 
central  core  of  this  fiber  can  be  fabricated  with  the  same 
material  and  properties  as  the  single-mode  fibers  (SMFs)  of 
the  laser  cavity,  easing  the  alignment  of  the  coupler  with 
the  cavity  fibers  and  reducing  splicing  losses  [16], 

The  fiber  profile  is  designed  such  that  it  supports  only  two 
zero-angular  momentum  guided  modes,  referred  to  as  the 
first-  and  second-order  guided  supermodes.  The  electric  field 
in  the  TCC  can  be  expressed  as  the  summation  of  these  guided 
modes  as 

E(x ,  y,  z,t)=  J2  Af  (s<  (•*•  y)ei(a^«*^ ,  (3) 

li= 1,2 


where  is  the  slowly  varying  envelope  of  the  electric 

field  of  the  pth  mode  with  the  normalized  spatial  distribution 
Fh(x,y)  and  propagation  constant  The  generalized  non¬ 
linear  Schrodinger  equation  describing  the  longitudinal  evolu¬ 
tion  of  a  pulse  in  this  fiber  coupler  can  be  written  as  [17,18] 


dz 


=  iS^A 


dt~l  2  df 


+  i 


(- ?)  E  .wvm;. 

'  6  '  v,K,p=l,2 


I*  =  1,2.  (4) 


The  indices  can  take  the  value  of  1  or  2,  corresponding  to  first- 
or  second-order  supermodes,  respectively,  pj}  *  denotes  the 
group  velocity  and  p\,  is  the  GVD  parameter  of  the  //th  mode. 
We  define  <5/?‘0)  =  =  pf-pf^  and  <5/?™  =  -Sp^  = 

p2  -  where  pf^  =  +  P ®)/2  for  i  =  0,  1.  The  nonlin¬ 

ear  index  coefficient  of  the  fiber  is  defined  by  n2  and  co0  is  the 
carrier  frequency.  The  nonlinear  coupling  coefficients  f^Kp 
are  given  by 


fpvKp  =  J  dxdyF;FvFKF (5) 

where  the  mutually  orthogonal  spatial  profiles  are  assumed  to 
be  normalized  according  to  ff  F\dxAy  =  ff  F\dxdy  =  1. 

In  the  two  different  TCC  fiber  couplers  studied  here,  in 
order  to  have  maximum  linear  coupling,  the  inner  and  outer 
radii  of  the  ring  core  are  designed  such  that  the  propagation 
constant  of  its  fundamental  mode  matches  that  of  the  central 
core.  The  radius  of  the  central  core  is  altered  from  its  original 
value  to  the  cutoff  value,  where  one  supermode  carries  70%  of 
the  input  power,  and  that  value  is  defined  as  the  maximum 
acceptable  fabrication  errors.  The  core  size  of  the  central  core 
is  the  only  geometrical  parameter  altered  in  this  paper,  be¬ 
cause  it  has  been  shown  that  the  performance  of  this  fiber 
is  more  sensitive  to  variations  of  the  central  core  radius  com¬ 
pared  with  that  of  the  ring  core  [19],  The  radius  of  the  central 
core  and  the  refractive  index  of  both  cores  in  coupler  A 
are  =  2.1  pm  and  v  t  =  1.465,  respectively,  whereas  the 
parameters  in  coupler  B  are  =4.1  pm  and  nt  =  1.451,  re¬ 
spectively.  The  operating  wavelength  is  2.0  =  1550  run  and  the 
cladding  refractive  index  for  both  couplers  is  assumed  to  be 
1.444  (corresponding  to  pure  silica).  Figure  2  shows  how  the 
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Separation  (pm) 

Fig.  2.  Maximum  acceptable  fabrication  errors  in  central  core  radius 
as  a  function  of  the  separation  between  the  cores. 

tolerance  to  fabrication  errors  in  these  two  couplers  varies  as 
a  function  of  the  separation  between  the  two  cores,  where  the 
separation  is  defined  as  the  gap  between  the  two  cores  (r2  -  rt 
in  Fig.  1_).  As  the  separation  between  the  ring  core  and  the 
central  core  increases,  the  overlap  of  the  fields  reduces,  which 
leads  to  a  smaller  coupling  and  a  larger  beat  length.  According 
to  the  results,  the  tolerance  to  fabrication  errors  reduces  as 
the  ring  core  is  placed  farther  away  from  the  central  core. 
Hence,  it  is  preferred  to  design  the  coupler  such  that  the  sep¬ 
aration  between  the  two  cores  is  small,  at  the  price  of  a  higher 
switching  power.  The  switching  power  is  defined  as  the  input 
power  at  which  the  output  power  from  the  central  and  ring 
core  are  equal.  Typical  fabrication  errors  in  core  radius  are 
around  ±  1%;  thus  according  to  Fig.  2,  the  maximum  allowable 
separation  between  the  two  cores  in  coupler  A  is  around  6  pm, 
whereas  in  coupler  B,  separations  up  to  10  pm  are  acceptable. 


whereas  high-intensity  noise  remains  in  the  launch  core,  re¬ 
sulting  in  a  high-intensity  contrast.  This  leads  to  the  formation 
of  optical  pulses.  The  power-dependent  transmission  feature 
of  the  coupler  can  be  used  for  mode-locking  fiber  lasers  [1,3]. 
In  modern  lasers,  the  saturable  absorber  is  a  critical  device 
only  for  allowing  the  pulse  to  build  up  from  noise  and  start 
the  lasing  process,  while  pulse  shaping  is  achieved  through 
using  a  narrow  spectral  filter  [21,22].  Therefore,  it  would  be 
reasonable  to  study  the  performance  of  the  device  for  nano¬ 
second  pulses,  which  are  pertinent  to  the  starting  process. 

A  hyperbolic  secant  optical  pulse  with  a  pulse  width  of  t0  is 
launched  into  the  nonlinear  coupler  through  the  input  SMF. 
The  nearly  Gaussian  spatial  mode  profile  of  the  input  beam 
excites  the  supermodes  of  the  coupler.  The  excitation  power 
related  to  each  supermode  can  be  found  using  the  overlap  in¬ 
tegral  of  that  particular  supermode  with  the  input  beam.  The 
pulse  propagation  inside  the  nonlinear  coupler  was  analyzed 
by  solving  Eq.  (4)  using  the  split-step  Fourier  method.  The  col¬ 
lected  output  power  from  the  coupler  is  Pout  =  \AC\2,  where  Ac 
is  the  field  amplitude  of  the  central  core,  obtained  by 

Ac(z,  t)  =  Ai(z.  t)  J F^Fydccdy  +  A2(z.  t)  jj F*F2dxdy ,  (6) 

where  Fi,  F2,  and  Fc  are  the  spatial  mode  profiles  of  the  two 
supermodes  of  the  TCC  coupler  and  that  of  the  propagating 
mode  of  the  central  core,  respectively.  The  finite-element 
method  was  used  to  calculate  these  profiles  [23].  Fc  is  calcu¬ 
lated  for  the  coupler  in  the  absence  of  the  outer  ring.  Trans¬ 
mittance  is  defined  as  the  energy  at  the  output  of  the  central 
core  divided  by  the  total  energy  injected  into  the  coupler, 
given  by 

T  =  f_+~\Ac(Lhb,t)\2dt 

/_+”(|A1(0,0|2  +  |A2(0,0|2)dt' 


3.  NONLINEAR  SWITCHING  BEHAVIOR 

The  TCC  fiber  can  be  used  as  an  all-fiber  mode-locking  device, 
referred  to  as  a  saturable  absorber,  as  shown  in  Fig.  3,  where 
the  beam  is  injected  into  the  central  core  of  the  TCC  coupler 
using  an  input  SMF  and  is  collected  at  the  output  from  the 
same  core,  using  another  SMF.  To  minimize  the  splicing  loss 
between  the  cavity  and  the  TCC,  the  central  core  is  considered 
identical  to  the  cavity  SMFs.  Given  that  the  injected  beam 
from  the  input  SMF  is  in  the  form  of  an  azimuthal  symmetric 
Gaussian  beam,  by  aligning  the  central  core  of  the  TCC  with 
the  input  SMF  core,  only  the  zero-angular  momentum  modes 
of  the  TCC  fiber  can  be  excited  [20],  which  is  the  case  con¬ 
sidered  in  this  paper. 

When  the  noise  of  the  gain  medium  in  the  laser  cavity  is 
transmitted  through  the  TCC  fiber  coupler  with  length  equal 
to  the  coupler’s  half-beat-length  Lhb,  low-intensity  noise  is 
coupled  linearly  to  the  ring  core  and  experiences  a  high  loss, 


nonlinear  switch 


input  SMF 

output  SMF 

>1 

m 

Fig.  3.  All-fiber  nonlinear  switching  device  composed  of  a  TCC 
coupler  placed  between  two  SMFs. 


In  Fig.  4,  the  normalized  input  pulse  and  output  pulses  cor¬ 
responding  to  different  peak  input  pulse  intensities  are  dis¬ 
played  for  coupler  A.  The  fiber  is  assumed  to  be  made  of 
silica  and  its  parameters,  with  reference  to  Fig.  1_,  are  rt  =  2.1, 
r2  =  8,  and  r3  =  9.95  pm.  These  parameters  are  chosen  such 
that  the  coupler’s  tolerance  to  fabrication  errors  is  around 
±  1%,  as  shown  in  Fig.  2.  It  can  be  observed  that  as  the  input 
power  increases,  a  larger  portion  of  the  pulse  energy  remains 
in  the  central  waveguide.  For  the  parameters  considered  here, 
a  switching  power  of  the  order  of  0.5  MW  is  required.  Tran¬ 
sition  from  low  transmission  to  high  transmission  occurs 
when  the  nonlinear  interaction  of  each  mode  with  itself  is 
comparable  to  the  linear  core-to-core  coupling;  therefore,  a 
fair  estimate  for  obtaining  the  power  threshold  for  nonlinear 
switching  in  an  ideal  coupler  without  any  fabrication  errors  is 
(n2w0/c)f  imPo  «  C,  where  C  =  <5/?0,i  is  the  linear  coupling 
coefficient  among  the  two  propagating  modes.  Hence,  a  lower 
switching  power  requires  a  coupler  with  a  large  nonlinear  co¬ 
efficient  and  a  small  linear  coupling  coefficient.  In  Fig.  5,  it  can 
be  seen  how  the  linear  coupling  coefficient,  C,  varies  as  a 
function  of  the  separation  between  the  cores  of  couplers  A 
and  B.  To  observe  better  the  difference  between  the  C  param¬ 
eter  of  these  two  couplers,  the  figure  is  plotted  on  logarithmic 
scale.  As  can  be  observed  in  the  figure,  for  small  separations, 
the  coupling  coefficient  in  coupler  A  is  larger  than  that  in 
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Fig.  4.  Normalized  pulse  power  at  the  output  of  the  central  core  for 
different  input  peak  powers  (units  are  MW)  when  pulse  width  is 
t0  =  1  ns.  A  similar  behavior  is  observed  for  picosecond  pulses, 
but  we  considered  nanosecond  pulses  because  they  are  more  relevant 
for  the  start  of  a  lasing  process. 

coupler  B,  whereas  vice  versa  occurs  for  larger  separations. 
Considering  a  fixed  separation,  the  effective  area  of  coupler  A 
is  less  than  that  of  coupler  B,  resulting  in  a  larger  nonlinear 
coefficient.  For  the  specific  couplers  considered  here,  the 
nonlinear  coefficient  of  coupler  A  is  approximately  3  times 
that  of  coupler  B  for  all  separations.  The  combination  of  linear 
coupling  and  nonlinearity  leads  to  a  lower  switching  power  in 
the  coupler  with  the  smaller  central  core  radius,  for  all  sepa¬ 
ration  values,  as  can  be  seen  in  Fig.  6.  To  reduce  the  switching 
power,  one  can  increase  the  separation  between  the  two  cores 
and  reduce  the  coupling.  However,  as  mentioned  in  the  pre¬ 
vious  section,  such  a  coupler  with  a  large  separation  between 
the  two  cores  requires  high  fabrication  precision  in  order  to 
operate  properly  in  the  linear  regime. 

In  Fig.  7,  the  sensitivity  of  the  coupler  transmission  has 
been  studied  against  ±  1%  variations  in  the  radius  of  the  cen¬ 
tral  core.  As  observed,  the  modulation  depth  of  the  nonlinear 
switch  is  reduced  in  both  cases.  When  there  are  no  variations 
in  the  central  core,  both  supermodes  are  excited  equally, 
whereas  in  STi  =  1%,  most  of  the  power  is  carried  by  the 
first-order  supermode  and  vice  versa  for  the  Sri  =  -1%  case. 
This  unequal  excitation  of  the  modes  is  the  reason  behind  the 
lower  quality  of  the  transmission  curves  of  these  two  cases 
compared  with  the  error-free  case.  Sri  =  1%  variations  have 


Fig.  5.  Coupling  among  the  two  cores  as  a  function  of  the  separation 
between  them  on  logarithmic  scale. 


Separation  (pm) 

Fig.  6.  Approximate  power  required  for  nonlinear  switching  in  cou¬ 
plers  A  and  B  as  a  function  of  the  separation  between  the  cores. 

the  advantage  of  shifting  the  transmission  curve  to  lower 
powers,  whereas,  5ri  =  -1%  variations  increase  the  power 
threshold  of  the  device,  which  is  not  desirable.  Larger 
amounts  of  error  reduce  the  modulation  depth  further;  thus, 
there  is  a  compromise  between  the  fabrication  error  and  the 
modulation  depth  of  the  device.  In  cases  where  a  lower  modu¬ 
lation  depth  is  sufficient  to  start  the  lasing  process  of  the  laser, 
one  can  deliberately  increase  the  central  core  radius  slightly 
from  the  designed  value,  and  hence,  reduce  the  switching 
power  required. 

4.  GROUP  VELOCITY 

The  two  supermodes  have  different  group  velocities;  thus,  the 
pulses  carried  by  these  two  modes  walk-off  due  to  traveling  at 
different  speeds  and  after  a  certain  propagation  distance, 
known  as  the  walk-off  length,  given  by  Lw  =  fo/(/?i,i  - 
the  two  pulses  will  not  have  any  overlap  and  hence,  will 
not  interact.  Given  that  the  length  of  the  coupler  is  chosen 
based  on  its  half-beat-length,  Lhb,  to  avoid  the  walk-off  effect, 
while  designing  a  fiber  coupler,  one  should  take  note  that  the 
minimum  pulse  width  that  can  propagate  through  the  coupler 
without  experiencing  walk-off  is  f0  =  Lhb  x  (fl\A  -  fii  z)-  For 
the  coupler  length  considered  in  this  paper,  the  minimum 
acceptable  pulse  width  that  will  not  experience  walk-off 
is  f0  «  14.7  fs. 


Fig.  7.  Relative  power  transmission  of  coupler  A  is  plotted  as  a  func¬ 
tion  of  input  power  with  variations  in  the  central  core  radius. 
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5.  GROUP  VELOCITY  DISPERSION 

When  the  two  cores  are  coupled  strongly  and  the  separation 
between  them  is  not  large,  the  GVD  and  also  dispersion  are 
comparable  to  that  of  SMF-28,  and  fluctuations  due  to  fabri¬ 
cation  imperfections  do  not  affect  the  dispersion  significantly. 
It  is  interesting  to  note  that,  although  for  large  separations  the 
dispersion  in  the  specific  fiber  geometry  considered  here  in¬ 
creases  rapidly  [24],  possible  fabrication  fluctuations  reduce 
the  dispersion  considerably.  Hence,  fabrication  uncertainties 
not  only  are  not  a  limiting  factor  in  terms  of  the  waveguide 
dispersion  but  in  fact  help  reduce  the  dispersion  as  well.  In 
addition,  for  the  pulse  widths  considered  here,  the  dispersion 
length  of  the  fiber  is  much  longer  than  the  coupling  length,  so 
even  high  amounts  of  dispersion  do  not  affect  the  coupler’s 
performance  considerably. 

6.  CONCLUSIONS 

From  a  practical  point  of  view,  we  have  studied  the  design 
considerations  for  multicore  fiber  couplers  used  in  nonlinear 
switching  and  mode-locking  applications.  The  effect  of  fabri¬ 
cation  imperfections  on  the  performance  of  a  two-core  fiber  is 
investigated,  which  can  be  generalized  to  all  multicore  cou¬ 
plers.  It  is  observed  that  there  is  a  compromise  between  fab¬ 
rication  error  and  the  required  nonlinear  switching  power, 
and  also,  the  modulation  depth  of  the  device. 

The  geometrical  parameters  considered  in  this  paper  are 
chosen  such  that  the  coupler’s  tolerance  to  fabrication  errors 
is  around  ±  1%.  This  limits  the  maximum  separation  between 
the  two  cores,  which  in  turn  increases  the  switching  power.  In 
a  more  precise  fabrication  process,  the  separation  between 
the  two  cores  can  be  increased  and,  consequently,  the  re¬ 
quired  switching  power  is  reduced.  For  example,  if  the  sepa¬ 
ration  is  increased  from  its  current  value  of  5.9  to  7.9  pm,  the 
switching  power  can  be  reduced  by  a  factor  of  3.5  down  to 
140  kW.  To  further  decrease  the  switching  power,  the  central 
core  radius  can  be  designed  such  that  it  is  slightly  larger  than 
its  optimum  value.  In  this  case,  in  order  to  achieve  the  re¬ 
quired  optical  properties,  the  errors  are  designed  directly  into 
the  fiber  preform.  This  lower  switching  power  is  achieved  at 
the  expense  of  a  lower  modulation  depth,  but  this  is  not  much 
of  an  issue  given  that  in  many  cases  even  a  low  modulation 
depth  is  sufficient  to  start  the  lasing  process  of  the  laser.  It 
should  also  be  pointed  out  that  to  decrease  the  power  thresh¬ 
old,  the  coupler  can  be  fabricated  from  material  with  a  higher 
nonlinearity  compared  with  silica  [25], 

Given  that  there  are  many  possible  ways  to  build  a  mode- 
locked  fiber  laser  cavity,  we  did  not  want  to  limit  ourselves  to 
one  particular  setup  and  so  we  did  not  discuss  any  specific 
laser  configurations  here.  Otherwise,  in  our  recent  publication 
[5],  we  have  shown  the  possibility  of  mode  locking  a  laser  cav¬ 
ity  with  multicore  optical  fibers.  Also,  a  mode-locked  fiber 
laser  using  a  waveguide  array  has  been  demonstrated  exper¬ 
imentally  [4],  which  serves  as  a  convincing  experimental 
proof  of  concept  for  the  potential  application  of  nonlinear 
mode  coupling  as  saturable  absorbers. 
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Abstract — A  detailed  analysis  of  four-wave  mixing  in  graded- 
index  multimode  fibers  is  presented.  We  derive  a  simplified 
analytic  expression  for  calculating  the  phase  mismatch  occurring 
as  a  result  of  waveguide  dispersion  among  the  pump,  signal  and 
idler  waves  in  these  fibers.  According  to  our  analysis,  we  have 
found  a  novel  expression  that  suggests  which  mode  combinations 
of  the  participating  beams  will  exhibit  identical  phase-matching 
conditions. 


I.  Introduction 

TO  overcome  the  capacity  limitation  of  nowadays  optical 
communication  networks  which  are  based  on  single  mode 
fibers  (SMFs)  [1],  [2],  many  researchers  are  looking  into 
using  multimode  fibers  (MMFs)  for  transmission  in  near  future 
networks.  The  reason  is  that  MMFs  can  further  increase  the 
per-fiber  capacity  through  employing  mode  division  multiplex¬ 
ing.  Consequently,  there  will  be  a  need  for  suitable  signal 
amplification  and  wavelength  conversion  schemes  in  these 
type  of  fibers.  Four-wave  mixing  (FWM)  in  these  fibers  is  a 
promising  technique  to  achieve  amplification  and  wavelength 
conversion  [3]— [8] . 

In  this  paper,  we  focus  on  FWM  in  graded-index  multimode 
fibers  (GIMFs).  Unlike  conventional  MMFs,  all  guided  modes 
in  a  GIMF  can  propagate  with  nearly  identical  group  veloci¬ 
ties  at  special  wavelengths  [9];  therefore,  nonlinear  coupling 
amongst  short  pulses  is  maximally  achieved  all  along  the  fiber, 
which  adds  to  the  richness  of  the  nonlinear  dynamics  [10], 
[11]. 

The  ultrashort  pulse  propagation  in  multimode  fibers  can 
be  described  by  a  generalized  nonlinear  Schrodinger  equation 
(GNLSE)  [12]— [14],  Using  the  GNLSE,  we  derive  a  set  of 
general  coupled  equations  describing  the  intermodal  FWM 
process  in  a  GIMF.  An  interesting  analytic  expression  regard¬ 
ing  the  phase  matching  condition  in  these  fibers  is  presented. 
For  different  spatial  mode  combinations  of  the  pump,  signal 


and  idler  the  phase-matching  conditions  are  analyzed  and  the 
nonlinear  coefficients  are  calculated. 


II.  Brief  overview  of  GIMFs 


In  this  section,  we  present  a  brief  overview  of  GIMFs  and 
establish  the  notation  that  will  be  used  in  the  rest  of  the 
paper  [15].  The  refractive  index  profile  of  a  GIMF  is  given 
by 


=  «o  (1  -  2A) 


P<R,  (1) 
P>  R, 


where  R  is  the  core  radius,  no  is  the  maximum  refractive  index 
in  the  center  of  the  core,  A  is  the  relative  index  difference 
between  the  core  and  cladding,  and  a  ~  2  characterizes  a 
near  parabolic-index  profile.  The  spatial  transverse  profile  of 
the  electric  field  of  a  propagating  mode  can  be  decomposed 
into  separate  radial  and  angular  components  as 

Fp,m(p,<jt)  =  Fp,m{p)eim<t>,  (2) 

where  [16] 

FpM  =  (3) 

P  o  zPo  P  0 

and 

JD 

n  —  11  A171  — 

P0  (4A2)!/4’  P 

L™  is  the  generalized  Laguerre  polynomial.  Here,  each  mode 
is  labeled  by  two  integers,  p  >  0  and  m  >  0,  referred  to  as  the 
radial  and  angular  mode  order,  respectively.  The  coefficients 
in  the  field  profile  in  Eq.  3  are  chosen  such  that  the  fields  are 
orthonormal  according  to 

p2tt  poo 

I  d(j)  I  pdpFpf  (p,  (/>)  —  (5) 

Jo  Jo 
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The  total  number  of  guided  modes  (counting  the  polarization 
degeneracy)  is  given  as 

Na  =  -^—nlk2R2  A,  (6) 

a  +  2 


where  A2  in  Eq.  3  is  also  given  by  Eq.  6,  and  k  =  2tt/X. 

The  near  parabolic  shape  of  the  refractive  index  profile  in 
Eq.  1  is  due  to  the  near  parabolic  Germanium  (Ge)  doping 
concentration  d  in  the  core  of  the  GIMF. 


d(p)  =  d0 
=  0 


P  <  A  (7) 
p>  R, 
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where  do  is  the  maximum  doping  in  weight  percentage  (wt%) 
at  the  center  of  the  core.  At  each  given  value  of  Ge  doping 
concentration,  the  refractive  index  can  be  accurately  described 
by  a  fit  to  the  Sellmeier  equation,  where  the  Sellmeier  coef¬ 
ficients  depend  on  the  Ge  doping  concentration.  In  particular, 
the  Sellmeier  coefficients  are  available  in  Refs.  [17],  [18] 
for  pure  fused  silica  (refractive  index  represented  as  nsi(A)), 
and  Ge-doped  fused  silica  with  a  Ge  doping  concentration 
of  6.3-wt%  (refractive  index  represented  as  716.3(A)).  Here, 
we  have  implicitly  assumed  that  the  refractive  index  of  the 
Ge-doped  fused  silica  increases  linearly  with  the  Ge  doping 
concentration  and  we  have  already  relied  on  this  assumption 
to  relate  Eqs.  1  and  7.  This  is  a  reasonable  approximation, 
which  is  also  supported  by  experimental  measurements  [19]. 
Therefore,  we  derive  a  generalized  Sellmeier  equation  to 
incorporate  the  dependence  of  the  refractive  index  on  both 
the  wavelength  and  Ge  doping  concentration  of  fused  silica  as 

n2(A,  d)  =  n2( A)  +  ^  [716.3(A)  -  n2(A)]  ,  (8) 


where  d  is  again  the  Ge  doping  concentration  in  wt%.  The 
utility  of  Eq.  8  is  that  it  provides,  along  with  Eqs.  1  and  7,  a 
reliable  estimate  of  the  refractive  index  for  a  broad  wavelength 
span  from  near-ultraviolet  through  mid-infrared  region  for 
common  levels  of  Ge  doping  concentration  in  GIMFs. 

All  the  modes  with  equal  group  mode  number  g  =  2p  + 

|  m  |  +  1  are  almost  degenerate  in  the  value  of  the  propagation 
constant,  which  is  given  according  to  the  formula 


where 


x  = 


2Ax)1/2, 

(9) 

2a/(a+2) 

(10) 

There  are  three  important  points  of  caution  regarding  the  use 
of  analytical  equations  presented  above.  The  first  point  is  that 
Eq.  3  is  only  strictly  valid  for  a  =  2;  however,  the  deviations 
are  small  and  Eq.  3  is  often  very  reliable.  The  second  point 
is  that  Eq.  3  is  only  a  solution  of  the  Maxwell’s  equations  in 
the  core  of  the  fiber  and  must  be  matched  to  an  exponentially 
decaying  modified  Bessel  function  of  the  second  kind  in  the 
cladding.  In  calculating  the  overlap  integrals  for  coupling  in 
and  out  of  the  GIMF,  as  long  as  the  beam-diameters  of  the  in¬ 
coupling  and  out-coupling  “often”-Gaussian  fields  are  smaller 
than  the  core  diameter  of  the  GIMF,  it  is  perfectly  fine  to  use 
Eq.  3.  For  larger  in-coupling  beam  diameters,  the  use  of  Eq.  3 
can  still  be  at  least  qualitatively  valid.  In  such  cases,  better 
accuracy  can  be  achieved  if  the  cladding  profile  is  properly 
accounted  for,  for  example  by  using  a  numerical  solver.  In 
either  case,  small  adjustments  in  the  normalization  coefficients 
of  Eq.  3  are  expected  according  to  Eq.  5  in  order  to  compensate 
for  substituting  Eq.  3  for  the  modified  Bessel  function  of  the 
second  kind  in  the  cladding.  The  third  point  is  that  the  presence 
of  the  cladding  removes  the  degeneracy  among  modes  within 
a  mode  group,  both  in  the  propagation  constant  and  the  group 
velocity.  While  this  effect  is  generally  quite  small,  it  is  more 
pronounced  for  higher  order  mode  groups,  which  have  a  larger 
overlap  with  the  cladding. 


III.  FWM  IN  A  MULTIMODE  OPTICAL  FIBER 

Four-wave  mixing  in  a  MMF  can  be  described  by  con¬ 
sidering  a  pump,  a  signal  and  an  idler  wave  oscillating  at 
frequencies  ojp,  ojs  and  jjt,  respectively.  During  the  FWM 
process  two  pump  photons  are  annihilated  and  two  photons 
at  the  signal  and  idler  frequencies  are  created  simultaneously, 
such  that  the  net  energy  and  momentum  are  conserved.  It 
is  assumed  that  all  three  waves  are  linearly  polarized  along 
the  x  axis  and  that  they  maintain  their  state  of  polarization 
while  propagating  in  the  MMF.  This  paper  focuses  on  FWM 
among  multiple  spatial  modes  and  therefore,  polarization- 
related  phenomena  are  left  for  future  work. 

The  (scalar)  electric  field  is  decomposed  into  a  sum  of  mode 
functions  Fp(x,y, oj)  with  envelopes  Ap(z,u>) 

E(x,y,z,t)  j  duJl~tN1^Y^A^z,uj)x  (11) 

M  v 

exp  [i(/3/j,(u})z  —  tit)]  +  c.c 

We  define 


(i//  (tu)  —  tf/i.  (A)  /3fun(Wref),  A  —  CU  CUref ,  (17) 

where  /3#i(cc)  is  the  propagation  constant  of  mode  y  at  fre¬ 
quency  oj,  and  /3fun(wref)  is  the  propagation  constant  of  the 
fundamental  mode  at  frequency  wref .  The  fundamental  mode  is 
taken  to  be  the  mode  with  the  largest  value  of  the  propagation 
constant  in  this  paper.  The  field  normalization  constant  Np(oj) 
is  chosen  such  that  the  mode  profile  satisfies  the  normalization 
and  orthogonality  condition 

JJ  dxdy  F*(x,y,u)F„(x,y,u)  =  6^.  (13) 

Considering  that  FWM  in  a  MMF  leads  to  a  wide  separation 
between  the  frequencies  of  the  pump,  signal,  and  idler,  it  is 
appropriate  to  use  separate  carrier  frequencies  for  each  of  these 
beams.  Therefore,  for  the  spatial  mode  y,  separate  envelopes 
Aj[,  A* ,  and  A^  are  assumed  for  the  pump,  signal,  and  idler, 
respectively.  The  upper  p,  s,  i  indices  refer  to  the  pump,  signal, 
idler  envelopes,  and  the  lower  indices  refer  to  the  spatial 
modes.  The  field  can  then  be  written  in  the  form  of 

Ap(z,t)  =  A^{zA)e^^)z  +  A*(2,f)einte^(^)z  (14) 
+  Alfl{z,t)e~mtei^uli'>z, 

where  fl  =  ojp  —  cos  =  uii  —  ojp.  It  is  assumed  that  <x>ref  =  ojp, 
which  leads  to  /3{un(ojp)  =  0.  The  optical  power  carried  by 
mode  y  is  \APi(z,t)\2,  where  A/i(z,f)  is  the  Fourier  transform 
of  Ap(z,oj). 

The  formulation  of  the  nonlinear  pulse  propagation  in  a 
multimode  optical  fiber  is  based  on  the  pioneering  work  of 
Kolesik  and  Moloney  [20]  and  Poletti  and  Horak  [12],  More 
recently  a  simplified  form  of  the  GNLSE  was  presented  in 
Ref.  [13],  which  describes  the  propagation  of  light  in  GIMFs. 
By  substituting  Eq.  14  in  the  GNLSE  of  Ref.  [13],  the  cou¬ 
pled  nonlinear  equations  governing  the  evolution  of  A(](z,t), 
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As  and  AlAz,t)  are  obtained  and  can  be  expressed  as 
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The  nonlinear  equation  for  the  idler  can  be  constructed  from 
Eq.  16  by  exchanging  all  the  signal  and  idler  indices  (s  o  i). 
These  equations  include  the  effects  of  self  phase  modulation 
(SPM),  cross  phase  modulation  (XPM)  and  FWM.  It  should 
be  noted  that  only  the  relevant  FWM  terms  that  satisfy  the 
conservation  of  energy,  2 up  =  ljs  ■  to,,,  are  retained.  In  the 
equations  above,  n2  is  the  nonlinear  index  coefficient  and  ‘Hz) 
are  the  spatial  phase  terms.  Three  representative  examples  are 
given  in  Eq.  17  from  which  all  the  other  phase  terms  can  be 
easily  written. 


*»)=exp 

®s™Kz(z)  =  exp 

<CM=exP 


H.-PZ  + PI +  &-%)* 
'iirPl  + %  +  &-%)* 


(17) 


pp  is  the  propagation  constant  of  the  spatial  mode  //  at  the 
pump  frequency  ujp:  fiP  =  and  so  on.  We  also  have 


J  =p,s,l. 


(18) 


Three  representative  examples  of  the  nonlinear  coupling 
terms  appearing  in  Eqs.  15  and  16  are  given  as 


llPZ  =  jj  dxdy  F;(cjp)Fv(up)Fk(u;p)F(*(up),  (19) 
1^%  =  jj  dxdy  F;(ua)Fv(ua)FK(ui)F^(Ui), 
lZ%  =  jj  dxdy  F;{wa)Fv(wp)FK((jp)FZ(Ui). 


all  other  nonlinear  coupling  terms  can  be  easily  written  using 
these  examples. 


IV.  FWM  IN  A  GRADED-INDEX  MULTIMODE  FIBER 

Equations  15  and  16  are  very  complex  and  include  the 
effects  of  SPM,  XPM,  and  FWM  processes.  In  general,  solv¬ 
ing  these  equations  requires  a  powerful  numerical  algorithm. 
However,  much  insight  will  be  obtained  by  making  a  few 
approximations  that  are  easily  justifiable  in  many  practical  sit¬ 
uations  [7].  In  the  following,  it  is  assumed  that  the  pump  field 
remains  much  stronger  than  the  signal  and  idler  and  remains 
undepleted  during  the  FWM  interaction.  This  simplifies  the 
equations  to: 


The  nonlinear  equation  for  the  idler  can  again  be  constructed 
from  Eq.  21  by  exchanging  all  the  signal  and  idler  indices 

(s  -o-  i). 

The  last  line  in  Eq.  21  contains  the  FWM  terms  that 
couple  the  pump,  signal,  and  idler.  Efficient  FWM  requires  the 
vanishing  of  the  effective  phase  mismatch  [7],  In  a  SMF,  the 
linear  phase  mismatch  is  often  compensated  by  a  nonlinear 
contribution  to  the  phase.  However,  in  a  MMF  such  as  the 
GIMF,  it  is  possible  to  satisfy  the  phase-matching  condition 
required  for  FWM  even  in  the  near-linear  regime  where  the 
pump  power  is  low.  The  linear  phase-matching  has  a  few 
distinctive  features  in  that  it  involves  multiple  modes,  signal 
and  idler  frequencies  are  far  from  the  pump,  and  the  effect 
of  the  nonlinear  phase  on  the  phase-matched  signal  and  idler 
frequencies  is  negligible.  The  FWM  phase-matching  condition 
requires 


(z)  =  1 


Afc  = 


p;+Pi-ppv-ppK  =  o 


(22) 


In  a  GIMF,  the  propagation  constants  of  the  modes  in  Eq.  9 
depend  on  the  mode  group  number.  Therefore,  the  phase¬ 
matching  condition  of  Eq.  22  can  be  cast  in  terms  of  the  mode 
group  numbers  of  the  pump,  signal,  and  idler  modes  expressed 
as  {ffp1^  ffp2'* ,  <7i,  gs}-  This  represents  phase-matching  between 
a  pump  beam  belonging  to  mode  group  qi! 1 ,  another  pump 
beam  belonging  to  mode  group  gp  ’ ,  the  idler  beam  belonging 
to  mode  group  g, ,  and  the  signal  beam  belonging  to  mode 
group  gs.  Eigenmodal  analysis  was  performed  using  the  finite 
element  method  in  order  to  calculate  the  propagation  constants 
of  the  beams  [21].  In  Fig.  1,  we  consider  possible  phase¬ 
matching  of  the  form  {1,1 ,  gi,  gs}',  therefore,  the  pump  modes 
are  assumed  to  be  in  the  lowest  order  mode  of  the  GIMF 
( gp 1')  =  gp2'*  =  1).  The  horizontal  axis  in  Fig.  1  represents  the 
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pump  wavelength  and  the  vertical  axis  represents  the  phase- 
matched  idler  wavelength.  The  solid  (black)  line  shows  the 
phase-matching  for  {1,1, 1,2},  where  the  two  pump  beams, 
idler,  and  signal,  propagate  in  the  modes  belonging  to  group 
numbers  1,  1,  1,  and  2,  respectively. 


pump  wavelength  (nm) 


pump  wavelength  (nm) 

Fig.  1.  The  phase-matching  condition  is  the  same  for  all  { I .  I  ,gi,gs } 
configurations  where  gi  +  gB  is  the  same.  Therefore,  {1,1,2,3}={1, 1,1,4}. 

V.  Phase-matching  in  graded-index  multimode 

FIBERS 

For  A  <  1,  Eq.  9  can  be  simplified  and  written  as 

Pg  ~  n0k  ^1  -  Ax  -  +0(A3),  (23) 

For  short  propagation  lengths  the  0( A2)  term  can  be  neglected 
because  it  does  not  amount  to  substantial  extra  phase  mismatch 
between  different  mode  groups.  For  situations  where  this  term 
can  be  ignored  and  a  «  2  can  be  trusted,  Eq.  23  transforms 
into 

V2A 

Pg  =  n0k - —  g,  (24) 

The  phase-matching  condition  is 

k  =  A  kM  +  Akw  +  A  kNL  =  0,  (25) 

where  k  is  the  summation  of  the  mismatch  occurring  due  to 
material  dispersion  (A kM),  waveguide  dispersion  (Akw )  and 
nonlinear  effects  (A kNL).  Assuming  that  A kNL  is  negligible, 
the  phase-matching  condition  is  satisfied  when  the  linear  phase 
mismatch  A k  =  A kM  +  A k^,  is  zero;  A kM  =  —A k^. 
The  mismatch  contribution  of  the  material  and  waveguide 
dispersion  for  {1,1, 1,3}  and  {1,1, 1,5}  configurations  are 


shown  in  Fig.  2.  Phase-matching  occurs  at  the  intersection  of 
the  A kM  and  —Ak^,  curves. 

Using  Eqs.  22,  24  and  25,  we  can  find 

A kM  =  n0(ujs)ks  +  n0(L0i)ki  -  2 n0(ujp)kp,  (26) 

Akw=-^jfiC™ M>  (27) 

where  is  defined  as 

Q™  =  9i  +  9  s  —  (ffp1)  +  ffp2^)  (28) 

Although  the  relative  index  step  A  is  a  function  of  the 
wavelength,  this  dependency  is  negligible.  Therefore,  A  can 
be  assumed  to  be  a  constant  value  for  all  wavelengths.  This 
approximation  is  considered  in  the  derivation  of  Eq.  27. 
Consequently,  according  to  Eq.  27,  AA ’w  can  be  approximated 
to  be  independent  of  the  wavelength  as  well  (seen  as  an  almost 
horizontal  line  in  Fig.  2),  and  so  it  only  depends  on  the  mode 
group  number  of  the  participating  waves.  Therefore,  different 
mode  combinations  of  pump,  signal  and  idler  that  exhibit  the 
same  C}VVM  value,  will  have  similar  phase-matching  conditions. 
This  means  that  if  in  separate  FWM  setups,  different  modes 
of  a  pump  beam  placed  at  a  fixed  wavelength  are  excited,  the 
FWM  interactions  which  have  identical  Q^^s  will  generate 
signals  (idlers)  with  identical  wavelengths. 


A(nm) 

Fig.  2.  The  phase-matching  curve  for  FWM  in  the  considered  GIMF  for  a 
pump  wavelength  of  850  nm.  The  A kM  curve  displays  the  contribution  from 
the  material  dispersion  and  the  A  kw  curves  display  the  waveguide  dispersion 
contributions  for  two  different  combinations  of  spatial  modes. 

For  a  better  understanding,  a  number  of  equivalent  configu¬ 
rations  are  mentioned  in  Table  1.  Although  all  of  these  spatial 
mode  combinations  satisfy  the  phase  matching  condition  and 
are  equivalent  wavelength-wise  (result  in  identical  signal  and 
idler  wavelengths  if  FWM  process  occurs),  they  have  different 
FWM  efficiencies.  The  reason  is  that  the  efficiency  depends 
on  the  overlap  integral  of  the  participating  spatial  modes  (see 
Eq.  19),  which  is  different  for  each  case.  This  can  be  observed 
in  Tables  II  and  III,  where  the  FWM  efficiency  of  some  phase 
matched  beams  are  displayed  through  their  nonlinear  coupling 
terms.  In  Table  II  the  pump  wavelength  is  Xp  =  850  nm,  while 
that  of  Table  III  is  1064  nm.  In  both  cases  it  is  assumed  that 
the  first  pump  beam  belongs  to  mode  group  =  1  which  is 
related  to  the  LPoi  mode.  Also,  it  is  worth  mentioning  that  the 
nonlinear  coupling  terms  will  only  have  non-zero  values 

when  —  to£  —  to^  +  to £  =  0,  where  is  the  angular 
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TABLE  I 

Mode  combinations  with  equivalent  phase-matching  conditions 


{1,1,1 ,2}={  1 ,2,2,2}={  1 ,2, 1 ,3} 
{1, 1,1, 3}={1, 1,2,2} 

{ 1, 1, 1,4}={  1, 1,2,3} 

{ 1, 1, 1,5}={  1, 1,3,3} 

(1,2, 1,2}={  1,3,2, 2}={  1,3, 1,3} 


mode  order  of  the  /rth  mode  at  frequency  j  (see  Eq.  2).  For 
instance,  the  combinations  {1,1, 1,2}  and  {1,1, 1,4}  and  their 
equivalents  do  not  satisfy  this  condition  and  thus  do  not  exist. 


TABLE  II 

A p  =  850  nm 


Pump  2 

Idler 

Signal 

rySppi 

LPoi 

LPoi  ,  683.9  nm 

LP02  ,  1122.8  nm 

0.608 

LPoi 

LP02  ,  684  nm 

LPoi  ,  1122.3  nm 

0.3645 

LPoi 

LPoi  ,  630.2  nm 

LP03  ,  1305.1  nm 

0.277 

LPoi 

LP03  ,  630.7  nm 

LPoi  ,  1302.9  nm 

0.0074 

LP02 

LPq2  5  683.7  nm 

LP02  ,  1123.1  nm 

0.2667 

LP2i 

LP21  ,  683.7  nm 

LP02  ,  1123.1  nm 

-0.0408 

LP2i 

LP02  9  683.7  nm 

LP21  ,  1123.1  nm 

-0.1649 

LPn 

LPn  ,  683.8  nm 

LP02  ,  1122.9  nm 

0.1346 

LPn 

LP02  9  683.9  nm 

LPn  ,  1122.7  nm 

-0.1075 

powers  at  the  end  of  a  fiber  of  length  L,  defined  as  Ps  and 
Pi,  are  expressed  as  [6] 

sinh2  ( gL ) 


Gs  =  ——  =  1  +  (7  sppiP0L)2 


(si)2 

a,  =  A  =  (Y™P,Lf  si,,h2(9i) 


P,. 


(gL)2 


(31) 

(32) 


where  Eqs.  31  and  32  describe  the  signal  amplification  and 
idler  generation,  respectively.  The  maximum  Gs  occurs  at 
the  signal  wavelength  where  the  phase-matching  condition  is 
perfectly  satisfied  (k  =  0).  As  k  increases,  g  decreases,  and 
drops  to  zero  for  k/2  =  ryspplP0.  If  k  is  further  increased,  g 
becomes  imaginary  and  for  gl  =  ±ii r,  Gs  =  0,  which  when 
substituted  in  Eq.  29  gives  us  the  bandwidth  as 


A  kb  =  2 


(33) 


VIE  Conclusions 

We  have  presented  a  detailed  theoretical  description  of 
FWM  in  GIMFs.  After  simplifying  the  analytical  term  de¬ 
scribing  phase  mismatch  in  these  fibers,  our  analysis  has 
resulted  in  a  useful  universal  expression  for  finding  separate 
mode  combinations  of  pump,  signal  and  idler  that  have  similar 
phase-matching  conditions. 


TABLE  III 
A p  =  1064  nm 


Pump  2 

Idler 

Signal 

rySppt 

LPoi 

LPoi  ,  733.9  nm 

LPo2  ,  1933.7  nm 

0.5225 

LPoi 

LP02  ,  739.2  nm 

LPoi  ,  1898.2  nm 

0.2014 

LP02 

LP02  ,  734.3  nm 

LP02  ,  1931.1  nm 

0.2225 

LP2i 

LP21  ,  734.3  nm 

LPo2  ,  1931.1  nm 

0.0439 

LP2i 

LP02  ,  734.3  nm 

LP21  ,  1931.1  nm 

-0.1217 

LPn 

LPn  ,  734.1  nm 

LP02  9  1932.4  nm 

0.1923 

LPn 

LP02  ,  736.8  nm 

LPn  ,  1913.7  nm 

-0.1229 

VI.  Parametric  gain  and  bandwidth 

In  order  to  use  FWM  in  MMFs  for  amplification  or  wave¬ 
length  conversion  applications,  it  is  required  to  know  the 
amount  of  gain  and  its  bandwidth.  To  understand  the  basic 
features  of  parametric  gain,  we  restrict  our  studies  to  the  case 
where  the  pump  beam  propagates  in  a  single  fiber  mode.  The 
parametric  gain  coefficient  g  is  defined  as  [6],  [7] 

g  =  V(YpplP0) 2  -  («/ 2)2,  (29) 

where  P0  is  the  incident  pump  power  at  z  =  0,  and  k  is  defined 
as 

k  =  Ak  +  (Ypps  +  lippi  -  27 pppp)Po ,  (30) 

For  simplicity,  we  consider  the  case  where  only  the  signal 
Pso  is  launched  into  the  fiber  together  with  the  pump,  and  the 
initial  idler  power  is  zero.  Using  the  coupled  wave  equations 
describing  FWM,  the  general  solution  for  the  signal  and  idler 
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